


\ T V7 



ii' *. 



CHOOL 

Cbi-.s-Booa 






NAVAL POSTGRADUATE SCHOOL 

Monterey, CalKornia 




THESIS 



PILOT STUDY ON THE APPLICABILITY OF 
VARIANCE REDUCTION TECHNIQUES TO THE 
SIMULATION OF A STOCHASTIC COMBAT MODEL 




by 




Curtis Smith 




September 1987 


Thesis 


Advisor: D. P. Gaver 



Approved for public release; distribution is unlimited 



T 234388 




SECURE q^Si/^ATlQN QF ThiS PAGE 



REPORT DOCUMENTATION PAGE 



i* REPORT SECURITY CLASSIFICATION 

UNCLASSIFIED 



lb RESTRICTIVE MARKINGS 



2 s SECURITY Classification authority 



2b DEClaSSiFiOTiON /DOWNGRADING SCh£OULE 



3 distribution/ availability of report 

Approved for public release; 

distribution is unlimited. 



4 PERFORMING 0RGAN7ATI0N REPORT NUM8£R(S) 



S MONITORING ORGANISATION REPORT NuM0£R(S) 



6 s NAME OF PERFORMING ORGANiSAT ON 

Naval Postgraduate Schoc 



6b OFFICE SYMBOL 
1 (If Spph^le) 



7s NAME OF MONITORING ORGANISATION 

Naval Postgraduate School 



U ADDRESS ( C# fy S tste snd ZiPCode) 

Monterey, California 93943-5000 



?b aD 0«E$S (C.fy. Utt* tndlirCodt) 

Monterey, California 93943-5000 



8* NAME OF FUNDING /SPONSORING 
ORGANISATION 



8b OFFICE SYMBOL 
(if spphcsbie) 



9 PROCUREMENT INSTRUMENT IOC N T»F iC A TiON NUMBER 



8c A0DR£SS(C»ry Sfjf* snd ZIP C ode) 



10 SOURCE OF FUNDING NUMBERS 



PROGRAM 
ELEMENT NO 



PRO/ECT 


Task 


WORK JNlT 


NO 


NO 


ACCESSION NO 



1 T.TiE (include Security Clswfic stion) Pilot Stu< _ _ 

Reduction Techniques to the Simulation of a Stochastic Combat Model 



P£«SONA>. AufHOHIS) 



SMITH, Curtis 



' 3j 'Y PE OF Rf PORT 

Master's Thesis 


’ 3b T-ME COVERED 


14 OATE OE REPORT (Ytti Month Oty) 


'S PACE COwNT 


FROM _ TQ 


1987 September 


109 



6 SlP^lC ventary notation 



COSATi CODES 


f ElD 


GROUP 


SuB GROUP 















18 SLI0JEC T T£ RMS (£or>tmuf 0 n reverie if neceissry snd identify by block number) 

Synchronization, Monotunicity , Antithetic 
Variates, Stratified Sampling, Crude 
Simulation, Simulation Efficiency, VRTs 



9 ABSTRACT (Continue on re*ene if neceusry snd identify by block number) 

This thesis investigates the applicability of VRTs to the simulation 
of stochastic combat models. Ways of measuring the efficiency of a VRT 
are explored. Antithetic variates and stratified sampling are applied 
to the simulation of a trivariate Markovian combat model. Means of 
programming the antithetic variates and stratified sampling to reduce 
the inherent variability of uncertainty in the output data of the model 
are illustrated. Response surface regression models are used to 
characterize the performance of the antithetic variates and stratified 
sampling in the Markovian combat model. 



D $ ‘ R'9U T<0N i AVAilABiliTy o* abstract 
CLnClASSiE'EOTJNL'MiTED □ SAME AS RPT □ QTi< oSt 



abstract security classification 

UNCLASSIFIED 



lit NA M£ ot RESPONSIBLE NOiViOUAC 

Prof. D. P. GAVER 



22b TELEPHONE (irxlud* ArttCodt) 

408-646-1605 



22 c OEE'CE SYMBOL 

55Gv 



00 FORM 1473. 8« MAR 



83 APR fd't'On msy b« ui«d until fiMwtttd 
All otn*f td t-O n\ sie ODiOlftt 



SECc»iTY C( ASSiFiCATiON HF Thi^ PAGE 



Approved for public release; distribution is unlimited. 



A Pilot Study 
Techniques 



on the Applicability of Variance Reduction 
to the Simulation of a Stochastic Combat 
Model 



by 

Curtis Smith 

Captain, United States Army 
B. S., United States Military Academy, 1 978 



Submitted in partial fulfillment of the 
requirements for the degree of 



MASTER OF SCIENCE IN OPERATIONS RESEARCH 



from the 



NAVAL POSTGRADUATE SCHOOL 
September 1 987 



ABSTRACT 



This thesis investigates the applicability of VRTs 
to the simulation of stochastic combat models. Ways of 
measuring the efficiency of a VRT are explored. 
Antithetic variates and stratified sampling are applied 
to the simulation of a trivariate Markovian combat 
model. Means of programming the antithetic variates and 
stratified sampling to reduce the inherent variability 
of uncertainty in the output data of the model are 
illustrated. Response surface regression models are 
used to characterize the performance of the antithetic 
variates and stratified sampling in the Markovian 
combat model . 
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INTRODUCTION 



I . 



A . GENERAL 

War policies and plans for military operations made 
during peacetime are significant for the mission 
accomplishments of combat operations conducted during 
wartime. Since experimentation with real combat is 
infeasible, military analysts use stochastic combat 
simulation models to study the effects of policy making 
on combat operations. The analysts’ inferences drawn 
from the results of these models are important to the 
decision maker since he has to use them to make the 
same decisions about military operations as he would if 
he could experiment with real combat itself. The output 
data from these models are realizations of random 
variables distributed around the values of the 
parameters of interest, or the models’ true 
characteristics, so the analysts can only estimate 
these parameters with error. The magnitude of error 
for each estimate can be measured in terms of precision 
or the variance of the estimate: if the estimate is 
unbiased then the smaller the variance, the greater the 
precision, and the smaller the error. Since the 
decision maker wishes to make decisions that are based 
on estimate(s) with a quantifiable error bound, the 
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analysts may find it possible, to apply specific 
statistical techniques to measure and control the 
variance of the estimate(s) to obtain a prescribed or 
at least quantifiable level of precision. 

The analyst’s capability of estimating a parameter 
of interest with high precision depends on the extent 
to which he is able to control the sample variance. 
When the analyst uses the mean value of the sample 
output data as the estimate of a parameter of interest 
and when individual samples are independent, the 
coefficient of the variance of the estimator is reduced 
by a factor of 1 /n where n is the sample size of the 
output data. A large sample size yields an estimate 
with a small variance and high precision. Multiple 
replications to obtain a large sample size in complex 
stochastic models can be prohibitively expensive in 
terms of resources like money, internal computer time, 
computer storage space, etc. This is especially true 
for large-scale, complex stochastic combat simulation 
models, which often require hours rather than minutes 
for a single computed replication. Since available 
computer time is a compelling constraint on military 
studies competing for scarce resources, the analyst is 
usually given an allocated amount of time to simulate 
his model. This specified amount of time may affect a 
desired level of precision of the estimate(s) that the 
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analyst wishes to obtain from the simulation. Since 
the analyst can only execute a fixed number of 
replications within this block of time, the sample size 
(number of replications) may not be large enough to 
achieve a variance small enough to give the analyst an 
acceptable precision for the estimate(s). Hence, the 
analyst must either accept the particular level of 
precision and error associated with such variance or 
apply other specific statistical techniques which are 
more likely to produce a smaller variance, and hence a 
level of precision with which he can feel more 
comfortable . 

An economizing scheme in simulation to reduce the 
variance of the estimator is to intentionally distort, 
control, and modify the random properties of the input 
variables in the simulation model. The output data 
resulting from the manipulation of these random numbers 
are random variables which are designed to be much 
closer together and more closely distributed around the 
true value of the model’s parameter of interest than is 
the case with simple random sampling. A sample 
distribution resulting from such a variance-reducing 
scheme has the same mean value but a potentially 
smaller variance than the distribution of the sample 
without the usage of this scheme. The different 
techniques for doing this scheme are called Variance 
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certain of 



Reduction Techniques (VRTs). The effects of 
these, when applied to a combat model, are the subject 
of this thesis. 



B. BACKGROUND 

VRTs were initially used to evaluate multi- 
dimensional integrals. They have since been applied to 
small Monte Carlo simulation problems but have not been 
extensively utilized in large complex stochastic 
simulation models. The utilization of these variance- 
reducing techniques in real-world combat simulation 
models is even less common. Consequently, limited 
examples of applications of VRTs in these simulation 
models are found in the literature. The major reason 
for this is because the performance of the VRTs is 
suspected to be uncertain and unpredictable. The 
analyst has no guarantee that the usage of VRTs will 
work all the time. Futhermore, he has no way to know 
beforehand how much variance reduction he will get from 
the application of VRTs whenever they are effective. 
However, VRTs, in our opinion, promise to be powerful 
and effective tools in simulation if the issues of 
their performance in specific simulations are 
understood. In this section we will describe the effect 
that they can have on simulation studies. In Chapter V 
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we illustrate their effectiveness for a particular 
combat model . 

The effectiveness of a VRT may be measured by the 
relative efficiency of the simulation in obtaining 
estimate(s) with the utilization of this scheme, to 
the efficiency of a simulation under the same 
conditions without the VRT. Efficiency as Handscomb 
(1969, p. 253) defines it is 

Efficiency = 1 / (Variance * Work). (1) 

Here "Work" generally refers to computing time. 
According to Handscomb, variance reduction succeeds if 
the VRT increases efficiency. From Equation 1, we see 
that a decrease in variance and/or work will increase 
efficiency. Hence, variance reduction in simulation is 
more than solely a decrease in the variance of the mean 
of the estimators. Handscomb( 1 969 ,p • 253) calls a 
technique variance-reducing if it "reduces the 
variance proportionately more than it increases the 
work involved" or "does not reduce the variance at all 
in the usual sense, provided that it saves work." The 
work involved in attaining estimates by simulation has 
many attributes. Hammersley and Handscomb ( 1 964 , p . 22) 
suggest that the number of simulation runs epitomizes 
this work. However, we can easily measure this same 



work in terms of computing cost or/and simulation time. 
For it is the availability of these factors that 
ultimately determine the precision of the estimators. 
Hence, an effective VRT may not only produce more 
precise estimates but also economize the time and costs 
associated with the simulation to obtain the level of 
precision for those estimates. The efficiency of VRTs 
will be discussed more fully in Chapter III. 

This potential saving in computer time has 
stimulated the interest of the United States Army 
Concept Analysis Agency ( CAA ) in the utilization of 
VRTs. CAA has studied the effectiveness of a VRT in two 
of its larger and more complex simulation models. The 
results of these studies were mixed (Johnson, Bates, 
and Graham, 1985). CAA therefore recommended the 
continuation of the studies to investigate the 
applicability of VRTs to reduce the inherent 
variability in large, complex, stochastic combat 
simulation models. 

C. PURPOSE AND OBJECTIVE 

The purpose of this thesis is to provide additional 
insight into the applicability of VRTs to stochastic 
combat models, and to provide a base for future studies 
in the application of VRTs to large-scale, real world, 
stochastic combat simulations. The objectives of this 



thesis are plainly to identify those VRTs that are 
applicable, and then to exhibit their performance in 
the applications to a class of stochastic combat 
simulation models. The question to be answered is: "Can 

VRTs be identified that are consistently effective for 
reducing simulation time and cost?" 

D. PROBLEM 

The problem for this thesis is to increase the 
efficiency of a stochastic combat simulation model 
utilizing VRTs in terms of (1) increased precision of 
the model’s estimates for an allocated amount of 
simulation time, and (2) reduced computer time for a 
predetermined level of precision. 

E. APPROACH 

In his doctoral dissertation, Andr^asson (1972) 
showed that variance reduction in queuing systems is 
influenced by (i) the transformation of random numbers, 
(ii) the structure and parameters of the simulation 
model, and (iii) the choice of the model response 
quantity. Condition (i) is an attribute of the VRTs. 
Conditions (ii) and (iii) are characteristics of the 
model. To solve the problem stated above, we 
investigate the effects of the parameters of a 



stochastic combat model, 



described in Chapter IV of 



this thesis, on variance reduction. We then use those 
results to formulate our approach to increase the 
efficiency of this model in terms stated in the problem 
above . 

F. ORGANIZATION OF THIS THESIS 

This thesis is organized into 6 chapters. Chapter I 
is the Introduction chapter. Chapter II reviews the 
literature of VRTs in simulation. Chapter III discusses 
ways of measuring the efficiency of a VRT and explores 
the tradeoffs of measuring for increased precision of 
estimation and reduced computer time. Chapter IV is 
concerned with the simulation of a stochastic combat 
model and the programming for variance reduction in the 
simulation model. Chapter V deals with the 
applicability and performance of VRTs in the simulation 
model. In Chapter VI, we make conclusions about the 
applicability of VRTs in stochastic combat models and 
provide recommendations about their use in larger and 
more complex, stochastic models that are used to study 



real-world combat systems. 



REVIEW OF LITERATURE 



II . 



A. INTRODUCTION 

The VRTs that we use to solve the problem stated 
in Chapter I of this thesis are antithetic variates and 
stratified sampling . But first we review the literature 
of variance reduction in simulation. This chapter 
concentrates on the practical applications of VRTs in 
simulation models. We present a brief summary of works 
of scholars and experts on this subject. We then 
describe the basic concepts of two VRTs that we feel 
are applicable to large-scale, complex, stochastic 
combat simulation models. It is these two VRTs whose 
performances we later exhibit in the combat model in 
this thesis. 



B. SUMMARY OF PREVIOUS WORKS 

In the last 15 years interest in VRTs in simulation 
has stimulated much activity on this topic in the 
Operations Research community. This section does not 
comprehensively review all works that have been written 
in the literature, but it presents a brief overview of 
the utilization of VRTs in simulation. The purpose of 
this section is to summarize some of the studies of the 
general applicability of VRTs in simulation. 
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Hammersley and Handscomb (-1964) reviewed many of 
the simplest ideas of variance reduction in simple 
Monte Carlo problems as they can be applied in the 
fields of Mathematics and Physics. Their most easily 
understood examples and outstanding successes were the 
evaluations of integrals and applications to particle 
physics. Handscomb (1969, p.252-262) later suggested 
that VRTs be adapted to simulation. He acknowledged 
difficulties in predicting the effectiveness of the 
techniques in particular situations, but he did 
propose, in practice, "... to proceed by more or less 
inspired trial and error, learning by experience which 
tools serve one best [or which techniques are 
effective]." He also stated that it may be much harder 
to tell how much variance reduction may occur in large 
and complicated simulation problems. These issues 
remain major concerns for one using VRTs in large, 
complex, stochastic simulation models. 

Moy (1969, pp . 263-288) adapted several VRTs to 
simulation and investigated their applicability to 
queuing systems. He concluded that VRTs were indeed 
capable of working in the simulation of queuing 
systems. Kleijnen (1974, Ch . Ill), who has written 
probably the most comprehensive and most referenced 
documentation on the subject of VRTs in simulation, 
discussed the relevant differences between sampling in 
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Monte Carlo problems and sampling in stochastic 
simulation models. He showed that VRTs may be adapted 
to accommodate these differences. Kleijnen also 
presented a detailed description and critical appraisal 
of six techniques so devised for utilization in the 
simulation of large complex systems. These VRTs are 
stratified sampling , importance sampling , selective 
sampling , control variates , antithetic variates , and 
common random numbers . These six sampling techniques 
have become the most well-known and popular VRTs in the 
literature . 

Other less-known VRTs, however, have been applied 
to simulation. McGraft and Irving (1974) survey some 18 
different techniques for implementation in large scale 
simulation problems. McGraft and Irving include a 
comprehensive listing of the characteristics, 
advantages and disadvantages, and criteria for 
applicability to large simulation models, and 
demonstrate the effectiveness of several of these 
techniques with a military simulation application. 

Many other articles and papers have been written on 
the subject of VRTs. There are too many of them to list 
in this thesis, but the survey ranges from specific 
techniques to more general methods in simulation 
experimentation. Some of the most recent papers written 
about the general applicability of VRTs in simulation 



are Nelson (1985) and Cheng (1986). Textbooks that 
illustrate the applications with simple but excellent 
examples of variance reduction in simulation are Gaver 
and Thompson (1973, Ch . 12), Fishman (1978, Ch . 3), Law 
and Kelton (1982, Ch . 11), Morgan (1984, Ch . 7), and 
Brat ley, Fox, and Schrage (1987, Ch . 2). 

C. DESCRIPTION OF VRTs USED IN THIS THESIS 

Moy (1969, p. 263-288) experimented with antithetic 
variates and stratified sampling and showed that they 
are indeed capable of signif icantly decreasing 
variability in the simulation of simple queuing 
problems. Likewise, we wish to achieve similar results 
when we apply them to the simulation of our stochastic 
combat model. We do this in Chapter V. In this section, 
we discuss the underlying conditions and fundamental 
concepts in the applicability of antithetic variates 
and stratified sampling in simulation. 

1 . Antithetic Variates 

The method of antithetic variates is relatively 
well-known in the literature of variance reduction in 
simulation (Kleijnen, 1974). It is one of the most 
useful VRTs because of its simplicity and general 
applicability. When the method of antithetic variates 
is used, the sampling process is modified by the 
manipulation of random numbers. A simulation run 



produces a response from the ' original sequence of 
random numbers {r}; then, a second simulation run 
produces an antithetic response from the sequence of 
the complementary random numbers { 1 -r > . The average of 
the two responses is an observation on the sample 
output data of the stochastic simulation model. The 
mean value of this sample is estimated as the parameter 
of interest. 

The variance of this estimate is reduced if the 

responses of the first and antithetic runs of each 

replication are negatively correlated. Besides the 

interchanging of the random numbers in each run, two 

other conditions must occur to produce negative 

correlation between the runs. First, each response must 

be a monotonic function of its respective random number 

stream; that is, large values in each stream of random 

numbers should have an opposite effect on the response 

than the small values, and vice versa . The second 

condition is that the responses to the events in the 

first run must be synchronized with the responses to 

the events in the antithetic run. Synchronization, 

defined by Kleijnen (1974, p. 193), occurs 

...if the i ’ th random number r^ generates [in the 
stream of the first run] a particular event (e.g., 
arrival of customer j) then in the antithetic run 
(1 - r^) should generate the same event (i.e., not 

the arrival of customer j ’ where j ’ ^ j and not a 
service time ) . 

18 



If the antithetic variate-s methodology, coupled 
with required conditions, is designed into the 
simulation, then the average of the two negatively 
correlated responses will tend to produce an estimate 
with a high degree of precision. That is, if by chance 
{r} yields a response above the value of the true 
parameter of interest, then { 1 -r } should yield a 
response below the value of the true parameter. When 
these responses are averaged, the deviations between 
the responses and the true parameter approximately 
offset each other resulting in relatively small net 
variability in the output data. This idea can be shown 



mathematically. Let 


*1 


be 


the response of the first 


run; X 2 


, the 


response 


of 


the 


antithetic run; and Y, the 


average 


of X! 


and X 2 . 












Y = (Xt 


+ X 2 


) / 


2 




VARy 


= 1 /4 * 


( VAR X1 


+ VAR X2 + 2 * C0V X1 jX2 ) 






- 1/2 


* ( 


1 + 


C0RR X1 , X2 ) * VAR X1 



Clearly, a negative CORR^i X2 reduces VARy. If 
C0RR X1 >X 2 equals, or is close to, -1 , then the VARy is 
mathematically zero or very close to it. Hence, the 
antithetic sampling is designed in simulation models so 
that the correlation between the pair of responses is 
as close to -1 as possible. 
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Monotonicity and synchronization must be 
designed into a simulation program for a particular 
model . Kleijnen (1974), Law and Kelton (1982), and 
Bratley, Fox, and Schrage (1987) are excellent 
references that discuss ways to do this. We discuss our 
design to achieve these two conditions for antithetic 
variates in our model in Chapter IV. As stated before, 
the method of antithetic variates is simple to 
implement and requires little to no extra computer 
time. Because of simplicity of this VRT , examples of 
its applications are illustrated in nearly every 
textbook that considers the subject of VRTs . 

2 . Stratified Sampling 

The stratified sampling technique, discussed in 
this section and applied to the simulation model in 
chapter V of this thesis, is a different version of the 
stratified sampling that Moy , Kleijnen and other 
experts on VRTs have adapted to simulation. Handscomb ( 
1969, p. 261) calls this particular version of 
stratified sampling another form of antithetic 
sampling. Andriasson (1972, p. 6) refers to it as an 
antithetic transformation. Gaver and Thompson (1973, 
pp . 585-586) name it stratification extending an 
antithetic idea. It is indeed stratification in that 
the sampling process is modified so that the range of 
random numbers is divided into two or more strata from 
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which the simulation runs produce- responses. It has the 
antithetic flavor in that the responses in all strata 
are averaged together to get an observation which is 
part of the sample output data. This technique is also 
similar to antithetic variates in that its estimator is 
an average of correlated responses (Gaver and Thompson 
1973, p. 586). Likewise, this estimator tends to have a 
smaller variance. 

In our review of this technique, we saw no 
necessary conditions, like those for the antithetic 
variates, for this technique to be successful in 
simulation. The design of stratified sampling into our 
simulation model in Chapter IV is similar to that one 
given in Gaver and Thompson (1973, p . 586 ) . 

D . SUMMARY 

An abundant amount of material has been written on 
the subject of variance reduction. Techniques used to 
reduce the variance in Monte Carlo problems have been 
adjusted to do the same in simulation models. The 
applications of VRTs in simulation have been 
illustrated in queuing systems and simple textbook 
problems but successful applications to larger, more 
complex, real-world stochastic simulation models have 
not been so amply reported. There is no guarantee that 
VRTs will work spectacularly for every situation in the 
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simulation, and when they do work it is necessary to 
estimate the magnitude of the variance reduction. Pilot 
tests are encouraged to help resolve these issues. 
Antithetic variates and modified versions of stratified 
sampling are two of the more simple and easily employed 
VRTs and will be applied to a stochastic combat model 
in Chapter V. 
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III. EFFICIENCY OF VARIANCE REDUCTION 



A. INTRODUCTION 

In the last chapter we reviewed some studies that 
involved VRTs. In this chapter we discuss the problem 
of measuring the efficiency of a VRT . Comparing 
variances of a parameter of interest obtained from the 
simulations with and without the use of a VRT 
respectively, on an ordinal scale, may reveal if the 
VRT works, but it provides little information about how 
well the VRT works. Clearly, a quantitative measure is 
more desirable. Therefore, the manner or scale on which 
the efficiency of a VRT is measured should provide as 
much information as possible on the performance of a 
VRT. In particular, it should provide at least some 
base to answering the following questions: 

(i) "Does the VRT work?" 

(ii) "If so, how great is the variance reduction in 
terms of increased precision for estimating 
the parameter of interest?" 

(iii) "How great is the variance reduction in terms 
of simulation time saved for estimating the 
parameter of interest?" 

(iv) "What are the tradeoffs, if any, between the 
potential increase in precision and the economy 
of simulation time when applying VRTs?" 

In the next section we examine two methods that are 
usually used in the literature to measure the 
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efficiency of a VRT . We evaluate them in terms of how 
well they answer the questions above. In the third 
section, we offer a third alternative which is a hybrid 
of the two previous methods for measuring the 
efficiency of a VRT. This third method, we think, 
answers all four questions above and is used to measure 
the efficiencies of the antithetic variates and the 
stratified sampling techniques whose performance is 
exhibited in this thesis. In the fourth section of this 
chapter we show how to use the third method of 
measuring the efficiency of a VRT to obtain the 
tradeoffs between increased precision and reduced 
simulation. The last section is a summary of this 
chapter . 



B. ASSESSMENT OF VARIANCE REDUCTION 



In 


the 


literature 


the 


efficiency of a 


VRT is 


usual ly 


measured by 


(1 ) 


a decrease in the variance 


( Method 


*1 ) 


or (2) 


the 


relative efficiency 


of a 


simulation 


to obtain 


an 


estimate using a VRT 


to the 



efficiency of the simulation using no VRT (Method #2). 
Henceforth, we refer to a simulation without the use of 
a VRT as crude simulation . 

Method #1 is well defined in Kleijnen ( 1974, pp . 
106-107). Kleijnen uses this method by defining the 
efficiency of a VRT as a percentage of reduction in varianc 
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Method #1 = (Var 0 - Var-| ) / Var 0 ) * 100% 



( 2 ) 



where Var 0 and Var i are variances obtained in the same 
amount of simulation time for crude simulation and 
simulation applying a VRT respectively. The measure of 
efficiency of a VRT which Kleijnen introduces may be 
interpreted as that portion of the variance which is 
not achieved by crude simulation but is obtainable in 
the same amount of simulation using a VRT. The sign of 
this portion determines whether the VRT increases or 
decreases the precision; a positive sign reveals an 
increase and a negative sign, a decrease. The magnitude 
of the portion indicates how much of the precision is 
increased or decreased respectively. With this method 
we can also see that the VRT has an identical effect on 
reducing simulation time for a prescribed level of 
precision as it does on increasing precision. Method #1 
provides answers to three of the questions stated in 
the last section, but it does not resolve the question 
of tradeoffs for increased precision and reduced time 
in a simulation using VRT. 

McGrath and Irving (1974, p. 295) use Method #2 to 
measure the efficiency of a VRT. They initially used 
this method, shown as Equation (3), to equate the 
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relative advantage gained in simulation time by using a 
VRT . 

Method M2 = 

Varg/Var-] * (Simulation Timeg )/( Simulation Time-]) (3) 

where subscripts 0 and 1 represent crude simulation and 
simulation applying a VRT respectively. The relative 
efficiency that McGrath and Irving used to measure the 
efficiency of a VRT results in a factor by which the 
efficiency of a simulation is increased or decreased by 
using a VRT. If the value of this factor is greater 
than one, then the VRT works; otherwise, it does not. 
The magnitude of this reduction is the actual value of 
the factor. For example, if the value of the factor is 
5, then the simulation applying the VRT can obtain an 
estimate in 1/5 the simulation time required by the 
crude simulation for the same precision level. Method 
M2 may be viewed either as the reduction in simulation 
time when both simulations are to achieve the same 
variance, or as an increase in precision when both 
simulations are run for the same amount of time. This 
method, like Method M 1, answers only the first three 
questions proposed in the first section of this 
chapter . 
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c. 



THE HYBRID METHOD 



Methods #1 and #2 measure increased precision at a 
fixed simulation time, or a reduced simulation time at 
a fixed level of precision. They do not, on the other 
hand, measure increased precision at a level of reduced 
simulation time, or vice versa; nor do they provide a 
means to explore such a possibility. In Chapter I we 
emphasize that variance reduction may increase 
precision and reduce simulation time. The efficiency of 
a VRT , in our opinion, should reflect both effects so 
that we can explore the tradeoff of any combination of 
precision and simulation time. Method #3 offers such 
possibility and answers all four questions in the 
introduction section of this chapter. It is a mixture 
of Methods #1 and # 2 . Method #3 has Kleijnen’s idea of 
reduction in variance and McGrath and Irving’s use of 
relative efficiency. We define the efficiency of a VRT 
as a relative efficiency (RE), as shown in Equation 4, 
and later define it in terms of increased precision 
(IP) and reduced time (RT). 



Method #3 = 

Efficiency-j / Efficiencyg (4) 

where Efficiencyg and Efficiency-) are the efficiencies 
of the crude simulation and simulation applying a VRT 
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respectively. The definition of the efficiency of a 
simulation, identified by Equation 1, is the inverse of 
the product of the sampling variance of the parameter 
estimate and the work. Henceforth, we equate work to 
simulation time, which is the total time of the 
simulation model to obtain a parameter of interest and 
a specified variance. Such time may be computed as n 
replications times T (average) simulation time per run. 
If k runs are in a replication, then simulation time 
equals the product of kn runs and T (average) 
simulation time per run. When these variables are 
substituted in Equation 1 , the efficiency equation 
becomes Equation 5: 



EFFICIENCY = 1 / ( Var * k 



n 



T) 



(5) 



If we are to increase the efficiency of a 
simulation using a VRT , then we must attempt to 
decrease one of the parameters in Equation 2. The 
variable k runs per replication is a constant of the 
VRT. Specifically, the antithetic variates constant k 
is two; stratified sampling constant k is three in our 
study (it can be more); and for no VRT, the constant 
value of k is one. The variable T is model dependent; 
that is, its value depends on the input parameters of 
the model. Attempts to decrease this variable may be 
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futile; futhermore Bratley, Fox, and Schrage (1987, p. 
48) point out that relatively little can be done to 
decrease it. We discussed the relationship between the 
variance (Var) and replications (n) in Chapter I. They 
are , in essence, the variables we wish to decrease. 
Throughout this thesis we interchange the phrases 
decrease in variance with increased precision and 
reduction in replications with reduced time 
(simulation). If we substitute Equation 5 into Equation 
4, we get Equation 6. Note since T (average) simulation 
time per run is the same for both simulations, it is 
left out of the equation. 

RE = Var 0 /Var -| * k 0 /k 1 * n 0 /n 1 (6) 

If the RE value in this equation is greater than 
one, then the VRT successfully increases the efficiency 
of the simulation and is said to be strong; otherwise, 
it is said to be weak. A strong VRT decreases the 
variance so that precision is increased and simulation 
time is reduced. A weak VRT, on the other hand, does 
not decrease the variance as well as a strong VRT; in 
fact, a very weak (or subversive) VRT may increase the 
variance, which causes a reduction in precision and 
necessitates an increase in simulation time. In most 
simulation models a VRT may be strong for certain 
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conditions and weak for other conditions. In this 
thesis, we look for such characterizations of the 
antithetic variates and the stratified sampling in the 
stochastic combat model we describe in the next 
chapter. In the next section, we define Equation(6) in 
terms of increased precision and reduced simulation 
time . 

D. TRADEOFFS OF GAINING PRECISION AND SAVING TIME 
Let us define IP as a decrease in variance, 



IP - (Var 0 - Var-L) / Var 0 



(7) 



and RT as a reduction in simulation time 



RT - (n 0 k 0 - n^) / n 0 k 0 . 



( 8 ) 



Then 



Var-L = (1.0 - IP) * Var 0 



(9) 



and 



n^k-L - (1.0 - RT ) x ngk 



( 10 ) 



If we substitute Equations (9) and (10) into 
Equation (6), we get Equation( 11 ) . 



RE = 1/(1. 0 - IP) * 1/(1. 0 - RT) 



(ID 
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Equation 11 defines the relative efficiency which 
we defines as Method #3 of measuring the efficiency of 
a VRT, in terms of increased precision and reduced 
time. This equation resolved the unanswered question 
identified as (iv) in the introduction section of this 
chapter. With this equation, we can examine any 
combination of IP and RT . For example, suppose we 
measure the efficiency of a VRT to have a RE value of 6 
for the same amount of simulation time (Hint: RT=0). 
Substituting these values into Equation (11), we get IP 
= 5/6 or 83.3% increased precision. 

Suppose we only need to increase the precision to 
75% instead of 83.3%, then we can substitute the values 
for IP=3/4 and RE=6 (RE should not change) into 
Equation 1 1 . We now get RT=1/3 (Note, we increase the 
precision 75% and reduce the simulation time 33.3%). 
Likewise, with RE=6 for the efficiency of the VRT, 
examples of other combinations are ( IP=2/3 , RT= 1 /2 ) ; 
(IP-1/2, RT-2/3); and ( IP-0 , RT=5/6 ) . In fact, we may get 
any combination of (IP,RT) between 0 and 5/6. Note, 
however, if we want to increase the precision beyond 
83.3% or 5/6, we will get an increase in simulation 
time. That is the tradeoff in terms of more increased 
precision. For example, we will have to increase the 
simulation time to 2/3 or 66.7% to accommodate an IP of 
90% for a RE value of 6. In short, the information 
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obtained from method #3 is that we can make an estimate 
more precise and save simulation time simultaneously. 



E . SUMMARY 

In the literature, there are generally two methods 
of measuring the efficiency of a VRT . Method #1 is a 
decrease in the variance; Method #2 is the relative 
efficiency of a simulation using VRT to crude 
simulation. Both methods may determine if VRT works in 
a simulation model. They also may indicate the 
magnitude of the variance reduction in terms of either 
increased precision for a fixed simulation time or 
reduced simulation time at a fixed level of precision. 
In this chapter, we introduced a third method of 
measuring the efficiency of a VRT. It is a hybrid 
between Method #1 and Method #2. Method #3 offers 
exploration into the tradeoffs of increasing precision 
and saving simulation time for any efficiency value of 
a VRT. 
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IV. SIMULATION OF A STOCHASTIC COMBAT MODEL 



A. INTRODUCTION 

In the last chapter we discussed how to measure the 
efficiency of a VRT to determine how much we may save 
in simulation time or/and how much we may increase the 
precision of a parameter obtained by crude simulation. 
In this chapter we show how we may apply VRTs to the 
simulation of a combat system. The combat model which 
we have chosen to simulate and to apply the VRTs is the 
BCD Markovian model developed in the doctoral 
dissertation of Abdul-Latif Rashid Al-Zayani (1986). A 
modified version of this model, formulated by 
Professor Donald P. Gaver , is in Appendix A. This 
stochastic model may seem very simple, but its 
simulation provides invaluable insights into the 
applicability of VRTs to stochastic combat model. 

Beside being stochastic, the BCD Markovian model is 
also discrete and dynamic in nature; hence, it is a 
discrete-event simulation model. We refer those readers 
who want to know about the nature of discrete-event 
simulations to simulation textbook such as 
Morgan( 1 984 ) , Law and Kel ton( 1 982 ) , or Bratley, Fox, 
and Schr age( 1 987 ) . In this thesis, we describe the 
simulation of the BCD model in terms of discrete 
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events. We describe, in detail, the crude simulation of 
the BCD model in the next section. In Section C, we 
show how we applied antithetic variates and stratified 
sampling to this simulation. We summarize the chapter 
in the last section. 

B. CRUDE SIMULATION OF THE BCD MARKOVIAN MODEL 

We discuss the crude simulation of the BCD model in 
four parts. First, we describe the combat scenario; 
second, we define the characteristics of the model; 
third, we explain the simulation of the combat process 
in the model; and finally, we discuss a FORTRAN 
simulation program written for the model. 

1 . The Combat Scenario 

As part of an air defense command, a wing of 
aircraft defenders is responsible for defending an area 
against a hostile air attack from a group of bombers. 
When detection of an incoming threat occurs a flight of 
D defenders is launched to engage B bombers making the 
attack. When the two groups are within aerial combat 
range the defenders seek a one-to-one combat engagement 
with the bombers at a rate i. Only one free defender 
can engage a free bomber in combat; a bomber will 
generally attempt to avoid any engagement with a 
defender. A combat engagement lasts until either the 
bomber is killed or the defender is killed. A defender 



34 



kills a bomber at a rate a, and a bomber kills a 
defender at a rate s. Hence, at any instant during the 
combat process, a defender is either free and searching 
to fight a bomber, fighting a bomber, or killed. A 
bomber is, likewise, either free and eschewing 
engagement with a defender, engaging in combat with a 
defender, or killed. The combat process is continued 
until either force is completely killed off or the 
duration of the combat period is terminated after T 
units of time. 

2 . The Characteristics of the BCD Model 

Hartman (1985, p. 2-18) characterizes the 
structure of a combat simulation model as combat 
entities , attributes , and events . We use these 
characteristics to simulate the combat process of the 
BCD model in the next subsection. Combat entities in 
the BCD model are free bombers, free defenders, and 
combat engagements. Each entity has attributes that 
describe a combat scenario. For the bombers, the 
attributes are the number of bombers and the rate a 
bomber shoots down a defender; for the defenders, the 
number of defenders and the rate a defender shoots down 
a bomber; and for the combats, the number of combat 
engagements and the rate that a bomber and a defender 
engage in combat. 
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Law and Kelton (1982, p. 4) define an event in 
a discrete-event simulation as an occurrence which 
changes the state of the system. The BCD model has five 
events. The first event is the initialization of the 
air battle. The initialization event governs the 
initial battle conditions. The next three events, are 
the interim events in the combat process. These events 
are (1 ) a combat between a bomber and a defender, (2) a 
defender killing a bomber, and (3) a bomber killing a 
defender. The occurrence of an interim event changes 
the state of the combat process at time t. The state of 
the combat process of the BCD model is represented by 
the trivariate-Markov process { B( t ) , C ( t ) , D( t ) ; t > 0 > ; 
where, B(t) is the number of free bombers at time t, 
C(t) is the number of combat engagements at time t, and 
D(t) is the number of free defenders at time t. As a 
Markov process, the combat process moves from state to 
state according to one-step transit probabilities that 
depend only on the current state. The fifth and last 
event in the combat process is the termination of the 
air battle. The termination event manifests the "end of 
the battle" conditions. The values of the state of the 
system at the occurrence of the termination event 
reflects the battle outcome. These values are the 
numbers of bombers and defenders that are alive at the 
end of this air battle. We will consider these values 
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as the parameters of interest in. the simulation of the 
BCD model . 

3 . Simulation of the Combat Process 

We simulated the combat process by maintaining 
a "bookkeeping account" of the changes in the state of 
the combat process as the events occur . The process 
begins in the initial state { B( t ) , C( t ) ,D( t ) ; t = 0 } with 
the initialization event being the commencement of the 
air battle. Henceforth, we let a value of B(t) equal b, 
a value of C(t) be c, and of D(t) be d. The interim 
events change the value of the state of the combat 
process as following: 



EVENT 

New combat 

Bomber kills Defender 
Defender kills Bomber 



STATE 

b-1 ,c+1 , d- 1 
b , c- 1 , d+1 
b+1 , c- 1 , d 



The combat process spends X(b,c,d) units of 
sojourn time in state (b,c,d) until another event 
occurs. The sojourn time X(b,c,d) is a random variable 
distributed exponentially with mean 

P(b,c,d) = 1 / ( iBD + (a + s)C) 
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for the time to the next event in- the combat process of 
the air battle, given that at time t the state was 
(b,c,d). Equation 12 is this sojourn time. To derive 
this expression, we use the inverse transform method to 
obtain unit-mean exponential random variables, where Vj 
is the j th random number in the sequence of a stream of 
uniform random numbers. The inverse transform method is 
discussed in the simulation textbooks listed in the 
reference section of this thesis. 

X(b,c,d) = - P(b , c , d) * ln(V-j) (12) 

We use the value of Equation 12 to advance the 
simulated time of the air battle as indicated by 
Equation 1 3 . 

t = t + X(b , c , d) (13) 

The combat process moves to another state when 
another event occurs. The probability of a specific 
interim event occurring is governed by an embedded 
Markov chain whose transition probabilities are as 
foil ows : 
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EVENT 



PROBABILITY 



New Combat 

Defender kills Bomber 
Bomber kills Defender 



fBD * P( b , c , d ) 
aC * P( b , c , d ) 
eC * P( b , c , d ) 



We again use the inverse transform method to 
obtain the conditions for an interim event to occur and 
to induce the change in the state of the combat 
process. Vj is the jth random number in the sequence of 
a different stream of random numbers. These conditions, 
events, and changes in the state of the combat process 
are listed below. 



CONDITION EVENT 

V j siBD* P( b , c , d ) New Combat 

V-; >iBD* P( b , c , d ) Defender kills Bomber 
and 

Vj < ( iBD+aC )#P(b,c,d) 

otherwise Bomber kills Defender 



STATE 

b - 1 , c + 1 , d+ 1 
b , c- 1 , d+1 

b+1 , c-1 , d 



Thus, we (i) generate a uniform random number 
to choose which interim event has occurred, (ii) update 
the state of the combat process, (iii) generate another 
uniform random number and transform it to an 
exponential random variable to determine the unit of 
time until the occurrence of the next event and to 
advance the simulated time of the combat process. We 
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repeat this procedure until the occurrence of the 
termination event. The termination event occurs when 
(1) all the bombers are killed (B(t) + C(t) = 0) or all 
the defenders are killed (D(t) + C(t) = 0) or (2) the 

time duration of the air battle has expired (t i T 

units of time). At the end of the aerial battle, the 

combat process is in state (b,c,d) from which we can 
compute the number of live bombers B (B(t) + C(t)) and 

the number of live defenders D (D(t) + C(t)). These are 

values of random variables for one possible battle 
outcome . 

4 . The FORTRAN Simulation Program 

We coded the crude simulation of the BCD model 
in FORTRAN. This FORTRAN program, consisting of a main 
program and four subroutines, is in Appendix B. The 
main program begins in an interactive mode. The program 
reads the values for the attributes of a combat 
scenario from the terminal and sends them to the BATTLE 
subroutine. BATTLE runs N replications of the combat 
process and returns the summary statistics of the 
outcome of N battles to the main program. The main 
program sends them to the STAT subroutine. STAT 
analyzes these battle statistics in terms of parameter 
estimates and returns the values of these parameter 
estimates to the main program. The main program then 
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sends these values of parameter estimates to a 
formatted output file. 

The BATTLE subroutine calls two subroutines 
UNIFOR and EXPON for the generation of U(0,1) random 
numbers. These two subroutines implement the 
congruential pseudo-random number generator 

u i+1 = 168071^ mod (2**31 - 1) (14) 

discussed and tested by Lewis and Orav (1985, Ch . V). 
UNIFOR generates a sequence of uniform random numbers 
for the selection of the occurrence of an interim 
event. EXPON generates a sequence of uniform random 
numbers for the computation of the unit of sojourn time 
in a state. 

The STAT subroutine performs statistical output 
analysis for the simulation. It computes the means and 
variances of the sample distributions of live bombers 
and defenders. The sample means for bombers and 
defenders 



B = SUM Bj_ / N (15) 

D = SUM Dj_ / N (16) 

are unbiased (point) estimators of E[B(t)] and E[D(t)] 
respectively. Similarly, the variances 
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(IV) 

(18) 



VARg = SUM ( B ± - B)**2 / N * (N-1) 

VARfj = SUM ( Di - D)**2 / N x (N-1) 

are unbiased estimators of VAR [E [ B( t ) ] ] and 
VAR[E[D(t)]] respectively (Larson, 1982). 

C. PROGRAMMING FOR VARIANCE REDUCTION 

In Chapter I, we noted that VRTs modify the 
sampling of random numbers. In this section, we discuss 
these modifications for the antithetic variates and 
stratified sampling in the simulation of the combat 
process. We describe the changes we made to the crude 
FORTRAN simulation program for the simulations using 
antithetic variates in Section 1 and stratified 
sampling in Section 2 respectively. 

1 . Antithetic Variates 

We make changes to the subroutines BATTLE, 
UNIFOR, and EXPON of the crude simulation program to 
use the antithetic variates. The FORTRAN program for 
the BCD simulation model applying antithetic variates 
is in Appendix C. The BATTLE subroutine computes the 
values of the parameters for one replication as the 
average values of the battle outcomes from a pair of 
runs of the combat process. We obtain the values of the 
battle outcome from the first run by using a stream of 
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uniform random numbers (U) and those of the second run 
by using a stream of complementary random numbers (1- 
U). Since we are attempting to decrease the variance of 
the estimates by inducing negative correlation between 
these two runs, we want to minimize this negative 
correlation. We, first, induce a negative correlation 
between the two runs by creating monotonicity between 
the random numbers and the values of the battle outcome 
within each run. We then minimize this negative 
correlation by synchronizing the sequences of random 
numbers (U) and the complement { 1 -U } (Bratley, Fox, 
Schrage 1987, p.47). Kleijnen (1974, p . 1 87 ) shows that 
a random variable generated by the inverse transform 
approach is monotonic. Hence we have monotonicity in 
the simulation since we used the inverse transform 
method to generate the uniform random variables in the 
simulation of the BCD model. 

Law and Kelton (1982, p. 352) indicate that the 
inverse transform approach also facilitates the 
maintenance of synchronization . With this method, we 
use only one uniform random number per sequence to 
obtain the desired random variable for each event in 
the combat process; as contrasted with other methods, 
like the rejection method, where we may use many random 
numbers to produce a single value for the desired 
random variable of the same events. Thus we initiate 
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synchronization in the model when we use the inverse 
transform method; nevertheless, we must still preserve 
it . 

We risk losing this synchronization in the BCD 
model because the number of interim events occurring in 
the combat process per run is a random variable. Hence 
the number of events occurring in the combat process in 
the first run may not be the same as the number of 
events occurring in the combat process in the 
antithetic run. Consequently, the number of random 
numbers needed in the antithetic run generally differs 
from that required in the first run. This phenomenon 
leads to the random number (U^) in the first run not 
being synchronized with the random number { 1 -Uj } of the 
antithetic run (Kleijnen 1974, p. 193). In other 
words, the complement of the j th uniform random number 
{ 1 -Uj ) is not used for the j th event in the combat 
process in the antithetic run. We are not able to 
control the random number of interim events in the 
combat process, but we can manage the way in which 
UNIFOR and EXPON generate uniform random numbers so 
that synchronization is maintained between the pair of 
runs per replication. 

We used the suggestions of Law and Kelton 
(1982, p.352) and Bratley, Fox, and Schrage (1987, 
p . 47 ) to maintain the synchronization that the inverse 
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transform method has initiated in the BCD simulation 
model. We modify subroutines UNIFOR AND EXPON to 
generate separate streams of random numbers {U} and 
complements { 1 -U } before simulating a pair of runs of 
the combat process. When the subroutine BATTLE calls 
UNIFOR and EXPON, it receives from each a two- 
dimensional array of random numbers, where the first 
column contains the stream {U} and the second column 
contains the stream { 1 -U } . Hence, we use the first 
column for the first run and the second column for the 
antithetic run. This approach guarantees that if the 
j th event in the first run uses (Uj), then the j th 
event in the antithetic run will use { 1 -Uj } . We do 
waste some of the random numbers in the arrays, but we 
do it judiciously. Since the number of random numbers 
used in the combat process is a random variable, we use 
only those random numbers that we need in each column 
and throw away the remaining so that no overlap is 
possible for the next pair of runs. As a result, we 
maintain synchronization. 

The last change we make to the crude simulation 
for the utilization of the antithetic variates is in 
the subroutine BATTLE. The subroutine BATTLE computes 
the values of parameters for each replication as 
follows : 
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Bi = + B 2 i) / 2 (19) 

where 

B 1 = the number of live bombers from the first run 

= the number of live bombers from the second run 

and 

Di = ( D 1 i + D 2 i) / 2 (20) 

where 

D - * = the number of live defenders from the first run 

D 2 ^ = the number of live defenders from the second run 

2 . Stratified Sampling 

As we stated in Chapter II, stratified sampling 
resembles the antithetic variates procedures, and so do 
the changes to the crude simulation. Hence we make 
changes similar to those in the simulation using 
antithetic variates for the simulation using stratified 
sampling. The FORTRAN program for the BCD simulation 
model using stratified sampling is in Appendix D. We 
modify subroutines UNIFOR and EXPON, where each 
generates a three-dimensional array of uniform random 
numbers from the three strata 



St = (0,1/3), S 2 = (1/3, 2/3), 



S 3 = (2/3,1) 



before simulating three runs of the combat process per 
replication. Note this does not need to be limited to 
3; we could have done more. 
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Gaver and Thompson (1973, p.586) describe this 
approach. First, a U(0,1) random number u-| -| is 
generated and placed in the first row and first column 
of the three-dimensional array. Next, 1/3 is added to 
the value of u-| i with a subtraction of one if needed to 
get u -|2 within the range of 0 and 1 . u-| 2 is placed in 
the first row and second column of the array. Next, 1/3 
is added to the value of u-| 2 » and if necessary 
subtracted by one, to get U 13 . u-j^ is placed in the 
first row of the third column of the array. If 
subroutine BATTLE calls for k random numbers, then k 
U(0,1 ) random numbers are generated, and the procedure 
obtains a value for each of the kx3 cells. BATTLE uses 
the first column of random numbers in the array for the 
first run, the second column for the second run, and 
the third column for the third run of the combat 
process . 

The values of the parameters for each 
replication are the average values of the battle 
outcomes from the three runs. The subroutine BATTLE 
computes the values of these parameters as follows: 



Bj. = ( B 1 i + B 2 i + B^) / 3 



(21 ) 



where 



B 1 ^ = the number of live bombers from the first r 



un 
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the number of live bombers from the second run 
the number of live bombers from the third run 



B 2 1 - 

B 5 i ■ 
and 

D ± = ( D 1 i + D 2 i + D5 i ) / 3 (22) 

where 

= the number of live defenders from the first run 

D 2 ^ = the number of live defenders from the second run 

= the number of live defenders from the third run 

D . SUMMARY 

The simulation of the BCD model is a discrete-event 
simulation. It begins with the initialization event and 
ends with termination event. The simulation of the 
combat process involves generating a sequence of U(0,1) 
random numbers to select interim event occurrences with 
changes in the state of the process and generating 
another sequence of U(0,1 ) random numbers to determine 
the unit of time until the next event occurs and to 
advance the simulated time of the combat process. 

The programming of the antithetic variates and 
stratified sampling modifies crude simulation. 
Monotonicity and synchronization are required in 
generating the uniform numbers for the simulation using 
these VRTs. Generating random numbers by the inverse 
transform method guarantees monotonicity. Generating 
sufficient random numbers by the inverse transform 
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method and in multi-dimensional arrays initiates and 
maintains synchronization. 
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V. APPLICATIONS OF THE ANTITHETIC VARIATES AND 

STRATIFIED SAMPLING 



A. INTRODUCTION 

We are now prepared to demonstrate the application 
of the variance reducing techniques to the simulation 
of a combat stochastic model. In this chapter we 
illustrate the performance of the antithetic variates 
(AV) and stratified sampling ( SS ) in the simulation of 
the BCD model. In Chapter IV we stated that the mean 
and variance of the parameters of interest estimated 
from simulation are used to analyze the output data of 
the model. Usually the estimated mean is of primary 
interest to decision makers, and the estimation of the 
variance is secondary. Since we use the variance of the 
parameters estimated from the simulation of the BCD 
model to exhibit variance reduction, we will, 
henceforth, focus on the variance. 

We examine the applicability of AV and SS in the 
BCD model by simulating many scenarios of the air 
battle and recording increases in simulation 
efficiency. We investigate AV and SS performance by 
mapping a response surface that characterizes the 
efficiency of variance reduction in the model. In the 
next section we specify the scenarios and discuss the 
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application of using AV and SS in the simulation of 
these scenarios. We then build response models that 
describe the performance of the AV and SS for these 
scenarios and discuss their experimental results in 
Section C. We present a brief summary of the chapter in 
the final section. 

B. APPLICABILITY OF ANTITHETIC VARIATES AND STRATIFIED 

SAMPLING TO THE BCD SIMULATION MODEL 

In Section B of the previous chapter we described 
the general scenario of the BCD model. In this section 
we specify various combat scenarios to observe how AV 
and SS are applicable to the simulation of the BCD 
model. Recall that seven attributes characterize a BCD 
scenario. Since a change of the values of one of these 
attribute will produce a different scenario, we chose 
to change the values for three attributes. We simulated 
10, 30 and 50 defenders against 10, 30, and 50 bombers 
at "end of battle" times of 25, 75, and 125. We 
maintain the a rate for a defender killing a bomber at 
.05, the 0 rate for a bomber killing a defender at .05, 
and the i rate for a bomber and a defender entering a 
combat engagement at .005. We also initialize every 
simulation with zero combat engagements. Thus, we 
observe the responses of the simulations at 3 "end of 
battle" times in 9 different scenarios. This 
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arrangement comprises a set of 27 independent 
simulations . 

We run this arrangement for crude simulation, 
simulation using AV , and simulation using SS. As a 
result, we perform a total of 81 different simulation 
experiments. In order to make a fair assessment of the 
applicability of AV and SS, we examine the variance 
obtained from the same amount of simulation or the same 
numbers of simulated battle runs for every simulation. 
We run 90 battles: this equates to 90 replications in 
crude simulation, 45 replications in simulation using 
AV , and 30 replications in simulation using SS. Table 
E.1 of Appendix E contains the statistical output for 
crude simulation; similarly, data in Tables E.2 and E.4 
are from simulations using AV and SS respectively. We 
use Equations 6 and 7 to measure the efficiency of the 
variance reduction (RE) and the increase in precision 
of the parameters from the simulation (IP) and place 
the AV results in Table E.3 and SS results in Table 



E.5. 



The values in Tables E.3 
SS respectively are applicable 
model . A RE value greater than 
VRT is effective in increasing 
positive IP value exhibits 
increase precision of the desir 



and E.5 show that AV and 
in the BCD simulation 
unity indicates that the 
simulation efficiency. A 
their effectiveness to 
ed parameter . With these 
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two values, we may also find the tradeoff of saving 
simulation time or gaining precision. We acknowledge 
that these values are all values of, or realizations 
of, random variables, but we believe that the tables 
show that the variance reduction adheres to a 
stochastic pattern. That is, the random variables 
obtained under certain scenarios will tend to have the 
same relationship to the random variables obtained 
under other scenarios. For example, the data in the 
tables suggest that high RE and IP values correspond 
to the scenarios that start combat with same numbers of 
bombers and defenders. The RE and IP values obtained 
under these scenarios appear consistently higher than 
the values under all other scenarios. Hence, the 
variance reduction measured by these RE and IP values 
are stochast ical ly greater than the variance reduction 
obtained from any other scenario. Since such even 
combats (i.e. equal combat power) are inherently more 
variable in outcome, the fact that variance reduction 
is greatest there is certainly welcome. In the next 
section we attempt to conduct a more thorough 
investigation of these phenomena so that we may 
understand how the AV and SS perform in the simulation 
of the BCD model . 
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C. PERFORMANCE OF ANTITHETIC VARIATES AND STRATIFIED 

SAMPLING IN THE BCD SIMULATION MODEL 

In the previous section we saw that AV and SS are 
applicable to the simulation of the BCD model. In this 
section we examine the variability of uncertainty in 
the model and then evaluate the applicability of AV and 
SS to reduce this uncertainty. We explore the changes 
in the AV and SS performance and examine the 
relationships of factors that affect these changes. 
Results of this analysis reveal the characterization of 
the AV and SS performance in the BCD model. 

1 . Experimental Design 

We use the data we generated in Appendix E to 
fit response surfaces that describe the uncertainty in 
the values of parameters in the BCD model and 

characterize the performance of AV and SS over the 
prescribed range of values in the three factors: 
initial numbers of Bombers and Defenders and "end of 
battle" Time. We code the three factors as 



x-j = (Time - 75) / 50, 

X 2 = (Defender - 30) / 20, 
x^ = (Bomber - 30) / 20. 



Each factor has 3 



levels . 



Thus, we may use a 3 5 
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factorial design to fit the data with the response 
surface equation 

E( y ) = eg + b-| Xi + + ei-|X-j**2 

+ e 22 x 2**^ "t * 2 + s i2^1*^2 

+ e-j^x-] *x-^ + S23 x 2 * x 3 (20) 

where 

E(y)= expected response 
60 = intercept 

6 ^= linear coefficient for factor i 
6 i i = quadratic coefficient for factor i 
e^j= interaction coefficient for the 
interaction of factors i and j 
x-^= level of factor i. 

We seek to obtain the maximum information from 
every observation; therefore, we chose a 3-^ Fractional 
factorial design. This design is the cuboctahedr on plus 
three center points (John, 1971). The three center 
points provide an unbiased estimate of error and 
repeated observations which permit us to test for Lack 
of Fit of the response surface equation we obtained. We 
use the cuboctahedron design to fit data for three 
response surfaces: ( 1 ) variability of uncertainty 

inherent in the battle outcomes, (2) efficiency of AV , 
and (3) efficiency of SS . This design, with its three 
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center points, and the data to which it is used to fit 
are shown in Tables 5.1 , 5.2, and 5.3. The variances 
data in Table 5.1 is the variability of uncertainty 
inherent in the battle outcome. We obtain this data 
from the variance of the estimate tabulated in the 
crude simulation table in Appendix E. The RE values 
data in Table 5.2 is the efficiency of AV for the 
estimation of the defender and bomber parameters. We 
took this data from the RE values in Table E.3. In 
Table 5.3 is the RE values data for the efficiency of 
SS . This data is obtained from the RE values in Table 
E.5. 



TABLE 5.1 DESIGN AND VARIANCES FOR A 3 x 3 EXPERIMENT 
ON THE VARIABILITY OF UNCERTAINTY INHERENT IN THE 
BATTLE OUTCOME 



DECODED LEVELS 


CODED 


LEVELS 


VARIANCES 


Time 


Defender 


Bomber 


x 1 


x2 


x3 


Defender 


Bomber 


25 


1 0 


30 


-1 


-1 


0 


. 0328 


. 0541 


1 25 


1 0 


30 


1 


-1 


0 


.0036 


. 1 725 


25 


50 


30 


-1 


1 


0 


. 1 687 


. 0826 


1 25 


50 


30 


1 


1 


0 


.3339 


.0359 


25 


30 


1 0 


-1 


0 


-1 


.0458 


. 0350 


1 25 


30 


1 0 


1 


0 


-1 


. 2247 


. 0045 


25 


30 


50 


-1 


0 


1 


.0628 


.1104 


1 25 


30 


50 


1 


0 


1 


. 0362 


.3479 


75 


1 0 


1 0 


0 


-1 


-1 


. 0387 


.0358 


75 


50 


1 0 


0 


1 


-1 


. 1 427 


. 0063 


75 


1 0 


50 


0 


-1 


1 


. 0064 


. 1 654 


75 


50 


50 


0 


1 


1 


.2717 


.2865 


75 


30 


30 


0 


0 


0 


.1618 


. 1895 


75 


30 


30 


0 


0 


0 


. 0925 


. 0966 


75 


30 


30 


0 


0 


0 


. 1 693 


. 1 589 
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TABLE 5.2 DESIGN AND RE VALUES FOR A3 x 3 EXPERIMENT 
ON THE EFFICIENCY OF THE ANTITHETIC VARIATES IN THE BCD 
MODEL 



DECODED LEVELS 


CODED 


LEVELS 


RE VALUES 


Time 


Defender 


Bomber 


x 1 


x2 


x3 


Defender 


Bomber 


25 


1 0 


30 


-1 


-1 


0 


2.5 


5.4 


1 25 


1 0 


30 


1 


-1 


0 


1 .5 


3.9 


25 


50 


30 


-1 


1 


0 


6 . 0 


2.2 


1 25 


50 


30 


1 


1 


0 


7.2 


1 . 4 


25 


30 


1 0 


-1 


0 


-1 


3.9 


3.0 


1 25 


30 


1 0 


1 


0 


-1 


4 . 0 


1 .5 


25 


30 


50 


-1 


0 


1 


2.8 


3 . 6 


1 25 


30 


50 


1 


0 


1 


1 .8 


10.9 


75 


1 0 


1 0 


0 


-1 


-1 


8. 1 


6 . 0 


75 


50 


1 0 


0 


1 


-1 


3.5 


1 . 1 


75 


1 0 


50 


0 


-1 


1 


1 . 1 


4 . 2 


75 


50 


50 


0 


1 


1 


10.1 


9.3 


75 


30 


30 


0 


0 


0 


6.4 


7.0 


75 


30 


30 


0 


0 


0 


13.5 


7 . 0 


75 


30 


30 


0 


0 


0 


10.8 


9.2 



TABLE 5.3 DESIGN AND 
ON THE EFFICIENCY 
MODEL 


RE VALUES FOR 
OF STRATIFIED 


A3 x 3 EXPERIMENT 
SAMPLING IN THE BCD 


DECODED 


LEVELS 


CODED 


LEVELS 


RE VALUES 


Time 


Defender Bomber 


x 1 


x2 


x3 


Defender 


Bomber 


25 


1 0 


30 


-1 


-1 


0 


2 . 6 


2 . 1 


125 


1 0 


30 


1 


-1 


0 


.9 


2 . 2 


25 


50 


30 


-1 


1 


0 


3.0 


1 . 9 


125 


50 


30 


1 


1 


0 


2 . 2 


1 . 2 


25 


30 


1 0 


-1 


0 


-1 


2 . 0 


2.8 


1 25 


30 


1 0 


1 


0 


-1 


2 . 6 


1 .7 


25 


30 


50 


-1 


0 


1 


2.9 


1 . 4 


1 25 


30 


50 


1 


0 


1 


1 .7 


2 . 6 


75 


1 0 


1 0 


0 


-1 


-1 


1 .6 


1 .5 


75 


50 


1 0 


0 


1 


-1 


1 . 6 


1 . 2 


75 


1 0 


50 


0 


-1 


1 


2 . 0 


1 .5 


75 


50 


50 


0 


1 


1 


3.0 


3.3 


75 


30 


30 


0 


0 


0 


3.7 


2 . 7 


75 


30 


30 


0 


0 


0 


2 . 9 


3.0 


75 


30 


30 


0 


0 


0 


2.5 


1 .5 
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2 . Experimental Analysis 



To perform the statistical analysis, we use the 
Response Surface Regression ( RSREG ) procedure in the 
Statistical Analysis System ( SAS ) computer software 
package on the IBM 370 mainframe. With this procedure 
we were able to obtain a second order response-surface 
equation by least-square regression, check for model 
adequacy, test for lack-of-fit, and identify critical 
surface values which were useful in helping to describe 
the shape of the surface. 

We fitted Equation 20 to the data in the 
respective design tables in the previous section and 
obtained multiple response surface equations and 
multiple analysis of variance (ANOVA) tables for 
corresponding responses. We assess the adequacy of each 
equation and test for fit from its corresponding ANOVA 
table. From each response surface equation, we 
generated additional data to obtain contour plots. We 
plotted contours of variability of uncertainty for the 
initial numbers of bombers and defenders at various 
Times. Similarly, we plotted contours of the 
efficiencies of AV and SS for initial numbers of 
bombers and defenders. From these plot we were able to 
see how the initial numbers of bombers and defenders in 
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combat affect the variability of uncertainty in the 
battle outcomes and the AV and SS performance. 

a. Variability of uncertainty of the battle 
outcomes 

Equation 21 provides an adequate 

description of the response surface that characterizes 
the variability of uncertain in the expected numbers of 
live defenders at the end of the battle in the BCD 
model . The response V;p e f' en cL er i s the expected amount of 
uncertainty in the defender estimate. 

= .1412 + .036xi + .104x2 — . 009 x 3 — . 014xi**2 

+ . 049 xiX 2 + . 008 x 2**2 - .05IX1X3 

+ .040x2X3 - .034x3**2 (21 ) 



TABLE 5.4 ANOVA FOR THE EXPECTED AMOUNT OF UNCERTAINTY 
IN THE DEFENDER ESTIMATE 



SOURCE 


d.f 


SS 


MS 


F-RATIO 


Fitted Surface 


9 


. 1 302 


.0145 


8.06 


Lack of Fit 


3 


.0102 


. 0034 


1 . 90 


Pure Error 


2 


.0036 


.0018 




Total 


1 4 


. 1 440 






R-Squar e = . 9043 


Mean 


Variance= .1194 


St d . 


Dev . = . 0525 



We see from ANOVA Table 5.4 that the 
variation in the variance of live defenders is 
insignificant at the 95^ level (its F-Ratio of 8.06 is 
less than F.g^.g 2 = 19.38). The Lack of fit is also 
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insignificant (F-Ratio= 1.90 < F .95;3,2 = 19.16). The 

R-Square value is .9043 which indicates that about 
of the variation in the variance of the expected number 
of live defenders is accounted for by Equation 21 . 
Further analysis reveals that the response surface is 
shaped like a rising ridge. The plots of contours at 
Figure 5.1 illustrate the nature of this Response 
surface. These pictures show that initial numbers of 
bombers and defenders affect the variance of defenders 
at Time 25. At Times 75 and 125 the initial numbers of 
bombers have little influence. Here the variance of the 
expected number of live defenders is affected solely by 
the increase in the number of initial defenders. Hence, 
as the number of defenders increases, the variability 
of uncertainty in the estimate of the expected number 

of live defenders at the end of the battle increases. 

Equation 22 provides an adequate 
description of the response surface that characterizes 
the variability of uncertainty in the expected numbers 
of live bombers at the end of the battle in the BCD 
model . V gombe r expected amount of uncertainty in 

the Bomber estimate. 

v Bomber = • 1 483 + .035 xt - .002x 2 + . 1 04x 3 - .031x 1 »*2 

- .041x-jX 2 + .032x 2 **2 - .067x-|X 3 

+ .038x 2 x 3 - .007x 3 **2 (22) 



60 



VARIANCE OF DEFENDER ESTIMATE AT TIME 25 




VARIANCE OF DEFENDER ESTIMATE AT TIME 75 




VARIANCE OF DEFENDER ESTIMATE AT TIME 125 




Figure 5.1 Contour Plots of the Response Surface for 
the Expected Amount of Uncertainty in the Defender 
Estimate 
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TABLE 5.5 ANOVA FOR THE EXPECTED AMOUNT OF UNCERTAINTY 
IN THE BOMBER ESTIMATE 



SOURCE 


d. f . 


SS 


MS 


F-RATIO 


Fitted Surface 


9 


.1331 


.0149 


6.77 


Lack of Fit 


3 


. 0073 


. 0024 


1 . 08 


Pure Error 


2 


. 0045 


. 0022 




Total 


1 4 


. 1 408 







R-Square= . 9 1 88 Mean Variances . 1 1 88 Std. Dev. =.0485 



ANOVA Table 5.5 shows that the variation in 
the variance of live bombers is insignificant at the 
95$ level (its F-Ratio of 6.77 is less than F. 95;9, 2 = 
19.38). The Lack of fit is also insignificant (F-Ratio= 
1.08 < F_95 j3 2 = 19.16). The R-Square value indicates 
that about 92$ of the variation in the variance of live 
bombers is accounted for by Equation 22. Further 
analysis reveals that this response surface is also 
shaped like a rising ridge. The plots of contours at 
Figure 5.2 illustrate the nature of this response 
surface. These figures show that the variance of the 
expected number of live bombers generally increase as 
the initial number of bombers increases, 
b. Antithetic Variates. 

Equation 23 provides an adequate 
description of an AV response surface in the estimation 
of the expected numbers of live defenders. 
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10 20 30 40 50 

INITIAL NUMGER OF BOMBERS 

VARIANCE OF BOMBER ESTIMATE AT TIME 75 




VARIANCE OF BOMBER ESTIMATE AT TIME 125 




Figure 5.2 Contour Plots of the Response Surface for 
the Expected Amount of Uncertainty in the Bomber 
Estimate 
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^Defender 3 10 -20 ~ -09x-| + 1.70X2 - . 46 x 3 - 4.25x-|**2 

+ .55x-|X 2 - 1.68x 2 **2 - . 28 x 1 X 3 
+ 3.40x 2 X3 - 2.85x3**2 (23) 

We see from ANOVA Table 5.6 that the variation in the 
efficiency of AV is insignificant at the 95% level (its 
F-Ratio of 1.46 is less than F^j.g 2 = 19.38). The 
Lack of fit is also insignificant (F-Ratio=.46 < 
F . 95 ; 3,2 = 19 . 16 ). The R-Square value is .8497 which 

indicates that about 85% of the variation in the 
Defender RE values is accounted for by Equation 23. 
Here, ^Defender is the expected simulation efficiency 
of AV generated to reduce the uncertainty in the 
Defender estimate. 



TABLE 5.6 ANOVA FOR THE EXPECTED EFFICIENCY OF AV (RE 
VALUE) GENERATED TO REDUCE THE UNCERTAINTY IN THE 
DEFENDER ESTIMATE 



SOURCE 


d.f . 


SS 


MS 


F-RATIO 


Fitted Surface 


9 


168.31 


18.70 


1 .46 


Lack of Fit 


3 


4.08 


1 .36 


.47 


Pure Error 


2 


25.69 


12.84 




Total 


1 4 


198.08 






R-Squar e= .8497 


Mean 


RE Value=5 . 54 


Std. 


Dev .=2.44 



Further analysis reveals that the response 
surface is shaped like a hill with a gentle slope on 
one side and a fairly steep slope on the other side. 
The maximum value of this surface occurs in the BCD 
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scenario that begins combat with 50 defenders and 40 
bombers and ends the battle at time 75. The plots of 
contours at Figure 5.3 illustrate the nature of this AV 
response surface. These pictures clearly show how AV 
performs in different scenarios for times of 25, 75, 
and 125. Beside having its best performance in a 
scenario that ends at time 75, AV appears to be strong 
in scenarios that initiate the air battle with at least 
30 defenders and 30 bombers. Its weakest performance 
seems to occur in those scenarios that commence combat 
with no more than 30 defenders and 30 bombers. The 
plots of contours show that the efficiency of AV is 
subversive in those scenarios whose simulation 
initializes the air battle with 30 or less defenders 
and 40 or more bombers. For these scenarios, simulation 
efficiency of AV may often increase, instead of 
decrease, the uncertainty in the Defender estimate. We 
will discuss why this is so in the Experimental Result 
Section. Similar analysis of the AV performance is made 
for the Bomber RE values. An adequate description of 
the AV response surface in the estimation of the 
bombers is characterized by 



Y Bomber = 7 *73 + .44 x-| - .69x 2 + 2 . 05x-j - 2.45x-|**2 

+ .18x-|X2 “ 2.05x 2 **2 + 2.20x-|X^ 

+ 2.50x 2 x-j - .53x3**2 
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(24) 




Figure 5.3 Contours Plots of the Responses Surface for 
the Expected Efficiency of AV Generated to Reduce the 
Uncertainty in the Defender Estimate 
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Y Bomber i s the expected simulation efficiency of AV to 
reduce the uncertainty in the Bomber estimate. The 
ANOVA Table 5.7 indicates that the proportion of the 
total variation in the Bomber RE values accounted for 
in Equation 24 is over 87 %. Furthermore, this variation 
is insignificant at the 95 % level (F-ratio value of 
8.19); lack of fit is also insignificant (F-Ratio= 
1.61). 



TABLE 5.7 ANOVA FOR THE EXPECTED EFFICIENCY OF AV (RE 
VALUE) GENERATED TO REDUCE THE UNCERTAINTY IN THE 
BOMBER ESTIMATE 



SOURCE 


d. f . 


SS 


MS 


F-RATIO 


Fitted Surface 


9 


118.74 


13.19 


8.19 


Lack of Fit 


3 


14.17 


4.72 


1.61 


Pure Error 


2 


3.23 


1.61 




Total 


1 4 


136.14 






R-Square= .8722 


Mean RE 


Value=5 . 05 


Std. 


Dev . = 1 .87 



Examining this response surface further we 
find that the shape of the surface changes over time. 
The plots of contours depicted in Figure 5.4 show that 
the shape of the surface looks like a saddle at Time 
25, a gentle slope at Time 75, and a uniformly rising 
ridge at Time 125. The critical values for this surface 
also change as its shape changes. Most notable are the 
values for maximum efficiency. At Time 25, maximum 
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EFFICIENCY OF AV IN SOMBER ESTIMATE AT TIME 25 





EFFICIENCY OF AV IN SOMBER ESTIMATE AT TIME 125 




Figure 5.4 Contour Plots of the Response Surface for 
the Expected Efficiency of AV Generated to Reduce the 
Uncertainty in the Bomber Estimate 
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efficiency occurs in scenarios that begin fighting with 
less than 20 defenders and bomber. By Time 125, maximum 
efficiency has shifted to the scenarios that start 
with at least 40 defenders and bombers. The least 
amount of AV reduction occurs, at any Time, in those 
scenarios that initialize the combat simulation with 
more than 30 defenders and less than 10 bombers. 

Here is a summary of what is revealed by 
the above analysis. AV , in general, seems to be the 
strongest and most consistent, and equally-distributed 
between closely-matched pairs of bombers and defenders . 
Furthermore, the larger the evenly-matched contest the 
greater the variance reduction. When the defenders and 
bombers are not evenly matched, AV is not as consistent 
and does not provide equal variance reduction in the 
estimation of the pair of parameters. It is strong in 
the estimation of the larger combatants and weak in the 
estimation of the smaller ones. 

c. Stratified Sampling. 

We analyze the efficiency of SS in the 
simulation of the BCD model in a similar manner as we 
analyzed the efficiency of AV . If we analyze 

Equations 25 and 26 in terms of ANOVA Tables 5.8 and 
5.9 respectively, we will get results similar to those 
we obtained in the last section. Therefore, we will 
forego this particular analysis. Note that Yj) e f en( ^ er is 
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the expected simulation efficiency of SS to reduce the 
uncertainty in the Defender estimate, and Y Bomt)er is 
the expected simulation efficiency of SS to decrease 
the uncertainty in the Bomber estimate. 

defender 3 3.03 “ -39 x 1 + .49x 2 + . 08x 3 - .15x 1 *»2 
+ .22x-)X 2 - .70x 2 **2 - .45 x-|X 3 + .55x 2 x 3 
- . 58 x 3 **2 ( 25 ) 



TABLE 5.8 ANOVA FOR THE EXPECTED EFFICIENCY OF SS (RE 
VALUE) GENERATED TO REDUCE THE UNCERTAINTY IN THE 
DEFENDER ESTIMATE 



SOURCE 


d. f . 


SS 


MS 


F-RATIO 


Fitted Surface 


9 


8.24 


.92 


2.49 


Lack of Fit 


3 


.53 


. 18 


.47 


Pure Error 


2 


.75 


.37 




Total 


1 4 


9.52 






R-Square= .86661 


Mean RE Value 


=2.27 Std. 


Dev .=.51 


^Bomber 3 2.90 - 


. 06x-j 


- . 0 3x 2 + . 


21 x 3 - .41xi* 


*2 


- .20X-|X 2 - . 


64x 2 **2 + . 


58 x 1 X 3 + .50x 


2 X 3 


- .56x3 


# # 2 






(26) 
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TABLE 5-9 ANOVA FOR THE EXPECTED EFFICIENCY OF SS (RE 
VALUE) GENERATED TO REDUCE THE UNCERTAINTY IN THE 
BOMBER ESTIMATE 



SOURCE 


d. f . 


SS 


MS 


F-RATIO 


Fitted Surface 


9 


5.18 


.58 


19.33 


Lack of Fit 


3 


1 .82 


.61 


20 . 33 


Pure Error 


2 


. 06 


.03 




Total 


1 4 


7 .06 







R-Square= . 7340 Mean RE Value=2.15 Std. Dev. =.61 

The contour plots at Figures 5.5 and 5.6 appear to 
have similar features. They show a relatively flat 
surface except at the corners. The corner with 50 
Bombers and 50 Defenders has the highest response and 
the other corners have low response. These plots 
suggest that the SS performance is generally consistent 
in all but a few scenarios in the BCD model. Maximum 
efficiency of SS occurs in those scenarios that 
initialize simulation with equally large numbers of 
bombers and defenders. It is very weak in those 
scenarios that begin combat with either less than 10 
defenders and more than 40 bombers or more than 40 
defenders and less than 10 bombers. 

5 . Experimental Results 

The experimental results can be summarized in 
Tables 5.10 and 5.11. Table 5.10 shows the 

relationships between the AV performance and the 
uncertainty in the Defender estimate and the SS 
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EFFICIENCY OF SS IN DEFENDER' ESTIMATE AT TIME 25 




EFFICIENCY OF S3 IN DEFENDER ESTIMATE AT TIME 75 




10 20 00 40 50 

INITIAL NUM3ER OF BOMBERS 

EFFICIENCY OF SS IN DEFENDER ESTIMATE AT TIME 125 




INITIAL NUMBER OF BOMBERS 



Figure 5.5 Contour Plots of the Response Surface for 
the Expected Efficiency of SS Generated to Reduce the 
Uncertainty in the Defender estimate 
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1 




10 20 00 40 50 

INITIAL NUMBER OF BOMBERS 

EFFICIENCY CF SS IN BOMBER ESTIMATE AT TIME 75 




EFFICIENCY OF SS IN BOMBER ESTIMATE AT TIME 125 




Figure 5.6 Contour Plots of the Response Surface for 
the Expected Efficiency of SS Generated to Reduce the 
Uncertainty in the Bomber Estimate 
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performance and the uncertainty in the Defender 
estimate. Similarly, Table 5.11 shows the relationships 
between the AV performance and the uncertainty in the 
Bomber estimate and the SS performance and the 
uncertainty in the Bomber estimate. 

We further examine this relationship by 
analyzing the data that measure the uncertainty (crude 
variance) and appropriate variance reduction (RE 
Values) in Appendix E. After we applied a logarithmic 
transformation to the data, we regress the RE values 
on the crude variance data and observe a strong 
logarithmic linear relationship between uncertainty and 
variance reduction. 



TABLE 5.10 RELATIONSHIP BETWEEN UNCERTAINTY AND THE 
EFFICIENCY OF VARIANCE REDUCTION IN THE DEFENDER 
ESTIMATE 



INITIAL 




UNCERTAINTY 


VARIANCE 


REDUCTION 


Defenders 


Bombers 


Variance 


AV 


SS 


10 


10 


medium 


strong 


fair 


10 


30 


medium 


fair 


fair 


10 


50 


smal 1 


weak 


weak 


30 


10 


1 arge 


strong 


fair 


30 


30 


1 arge 


strong 


strong 


30 


50 


medium 


fair 


fair 


50 


10 


1 arge 


strong 


weak 


50 


30 


large 


strong 


fair 


50 


50 


1 arge 


strong 


fair 
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TABLE 5.11 RELATIONSHIP BETWEEN UNCERTAINTY AND THE 
EFFICIENCY OF VARIANCE REDUCTION IN THE BOMBER ESTIMATE 



INITIAL 




UNCERTAINTY 


VARIANCE 


REDUCTION 


Defenders 


Bombers 


Variance 


AV 


ss 


10 


10 


medium 


strong 


fair 


10 


30 


large 


strong 


fair 


10 


50 


large 


strong 


fair 


30 


10 


smal 1 


fair 


fair 


30 


30 


large 


strong 


fair 


30 


50 


large 


strong 


fair 


50 


10 


smal 1 


weak 


weak 


50 


30 


medium 


fair 


fair 


50 


50 


large 


strong 


fair 



This relationship is manifested in the 
multiplicative equation shown in Table 5.12. VAR is the 
value of uncertainty or variance of the corresponding 
estimate obtained from crude simulation, and Y is the 
simulation efficiency of the variance reduction or RE 
value for the corresponding estimate. 

TABLE 5.12 ANALYSIS OF THE RELATIONSHIP BETWEEN 
UNCERTAINTY AND EFFICIENCY OF VARIANCE REDUCTION IN 
PARAMETER ESTIMATES 



VRT 


ESTIMATE 


RELATIONSHIP 


CORRELATION 

COEFFICIENT 


SET 

POINT 


AV 


Def ender 


Y= 10.43 


M 


VAR** . 362 


.8609 


. 00154 


AV 


Bomber 


Y = 9.73 


H 


VAR** . 362 


.8217 


.00186 


SS 


Def ender 


Y = 3.18 


H 


VAR** .144 


. 5601 


. 00032 


SS 


Bomber 


Y= 2.99 


M 


VAR** . 147 


.6054 


. 00058 
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The correlation coefficient reveals the 
strength of the logarithmic linear relationships 
between uncertainty and the efficiencies of AV and SS. 
With the values of the exponent in the equations being 
less than one and the values of the correlation 
coefficient being positive, the efficiencies of AV and 
SS are observed to increase, at a decreasing rate, as 
the uncertainty (variance) increases. We obtained the 
Set Point by setting Y=1 in the corresponding equation 
and solving for VAR. At this value, simulation 
efficiency nether increases nor decreases. Now if we 
observe a value of uncertainty, or variance obtained 
from crude simulation, above this set point, then we 
expected to get an efficiency of a VRT to increase the 
simulation efficiency. On the other hand, if the value 
is below the set point, then we expect the efficiency 
of the VRT to decrease the simulation efficiency. 

Here is the bottom line on AV and SS 



performance 


1 . If 


we 


of 


the 


a . 


We 




in 




(1 



in the BCD model : 



apply antithetic variates to the simulation 
BCD model , 

may expect the variability of uncertainty 
the defender estimate 

) to decrease if the variance of the 

estimate obtained from crude simulation 
is at least .00154, and 



76 



(2) to decrease, at a decreasing rate, with 
an increase in the variance of the 
estimate obtained from the crude 
simulation . 

b. We may expect the variability of uncertainty 
in the bomber estimate 

(1 ) to decrease if the variance of the 
estimate obtained from the crude 
simulation is at least .00186, and 

(2) to decrease, at a decreasing rate, with 
an increase in the variance of the 
estimate obtained from the crude 
simulation . 

2. If we apply stratified sampling to the simulation 

of the BCD model , 

a. We may expect the variability of uncertainty 
in the defender estimate 

(1) to decrease if the variance of the 
estimate obtained from crude simulation 
is at least .00032, and 

(2) to decrease, at a decreasing rate, with 
an increase in the variance of the 
estimate obtained from the crude 
simulation . 

b. We may expect the variability of uncertainty 
in the bomber estimate 

(1 ) to decrease if the variance of the 
estimate obtained from the crude 
simulation is at least .00058, and 

(2) to decrease, at a decreasing rate, with 
an increase in the variance of the 
estimate obtained from the crude 
simulation . 

E . SUMMARY 

We illustrated the performance of the antithetic 
variates and stratified sampling in the simulation of 
the BCD model. The manifestation of their pair 
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performance was characterized by response surface 
equations and plots of contour lines. Both VRTs were 
shown to be effective in increasing simulation 
efficiency, but they perform differently in the BCD 
model. AV provides the largest amount of variance 
reduction but is more volatile. AV increases the 
simulation efficiency on the average of 5 times the 
crude simulation; it is strong in the BCD scenarios 
where there is large amount of uncertainty in the 
battle outcomes for live bombers and defenders, and 
weak in those scenarios where there is little amount of 
uncertainty in the battle outcomes. SS , on the hand, 
has a more consistent performance. SS increases the 
simulation efficiency at a mean of 2 times the crude 
simulation. It performs nearly the sarnie in every 
scenario except where the uncertainty is large. 
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VI . 



CONCLUSIONS AND RECOMMENDATIONS 



A. CONCLUSIONS 

The objective of this pilot study has been to 
investigate the applicability of VRTs to reduce the 
inherent variability in stochastic combat models. We 
examined the effects of applying AV and SS to the 
simulation of a simple stochastic combat model . We have 
now shown that AV and SS are applicable to this 
stochastic combat model. We can infer that these VRTs 
are indeed capable of working in stochastic combat 
models, and their prospects in larger and more complex 
stochastic combat models are even more promising. The 
conditions of monotonicity and synchronization are 
essential parts of the design of the simulation program 
for these models. Hence, we feel that sizable increase 
in simulation efficiency is possible if these 
requirements are met in the simulation. 

The experimental results of applying VRTs to the 
BCD simulation model show that the strength of AV and 
SS is influenced by uncertainty. A strong variance 
reduction results from a large variance of the estimate 
obtained from crude simulation. A weak variance 
reduction is caused from a small variance of the 
estimate obtained from crude simulation. Hence, sizable 
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and consistent variance reduction depends on large 
variability of the simulated output from the stochastic 
combat model. Therefore, the variability of the output 
from larger and more complex stochastic models must 
also be large enough to obtain the size and consistency 
of simulation efficiency and variance reduction one 
desires from the applications of these VRTs to such 
models . 

B. RECOMMENDATIONS 

The pilot study presented in this thesis provides 
a base for further studies in the applications of AV 
and SS to large-scale, real world, stochastic combat 
simulation models. Usually complex simulation models 
have many subroutines or modules. The variability of 
uncertainty in the output data from these modules may 
vary from low to high. We recommend that a study of 
this matter focus on the degree of variability of 
uncertainty in the output data from each module. The 
interest of the study should be concerned with the 
relationship between the performance of the VRT and the 
variability of the output data from each module. The 
results should indicate where and how the VRT may be 
used in the model in order to maximize the simulation 
efficiency of the model. For example, if the study 
shows that a VRT performs strongly in a particular 
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module whose variability of output data is large and it 
performs poorly in another module whose variability of 
output is small, then the study should recommend that 
the VRT be used in the module in which it performs best 
and not be used in that module for which it performs 
poorly. Using VRTs in the module with small variability 
would most likely decrease simulation efficiency for 
that model and, at worst, suboptimize the overall 
performance of the VRT in the complex combat model. 
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APPENDIX A 



FORMULATION OF THE BCD MARKOVIAN MODEL 



by 

Professor Donald P. Gaver 
Naval Postgraduate School 
Monterey, California 



B bombers are approaching a group of D defenders. 
When the two groups approach within range each defender 
searches for a bomber; after he finds one they engage 
in combat. Either bomber or defender may win the 
combat; the survivor becomes "free", and is a candidate 
for the next combat. In general, bombers attempt to 
avoid combat, defenders seek it out. 

This situation becomes a tri-variate Markov chain 
if the following state is defined: { B ( t ) , C ( t ) , D( t ) > . 

Here t is conveniently measured from the time bombers 
and defenders are close enough to permit combat at all, 
B(t) is the number of free bombers at time t 
thereafter; ditto for D(t), the number of free 
defenders; C(t) is the number of one-on-one combats. 
Here are a set of transition rates: 

(1) Combats begin . If B(t)=b, C(t)=c, D(t)=d and if ? 
is the rate at which free bombers are found by 
free defenders then with probability 



82 



1 - ibdt + o( t ) 



(( e -*t)b)d = 

no defender finds a bomber in time t. Hence the 
probability that a defender does find a bomber is 
ibdt + o(t). This is the rate at which free 

bombers and defenders get converted to combats: 
the state 

(b, c, d) •» (b-1 , c+1 , d-1 ) with prob ibdt . 

(2) Defenders win . Same initial conditions. If a 
bomber is in combat with a defender the 
probability that a defender shoots down the 
bomber in time t (combat duration) if the latter 
doesn’t hit the defender is 

P{Combat duration s t I bomber doesn’t hit) 

= 1 - e -at . 

Likewise, the probability that the bomber shoots 
down the defender is 

P{Combat duration i t ! defender doesn’t shoot) 

= 1 - e -9t . 



83 



Now suppose they both shoot* doing so 
independently. Model the probability that both 
survive to t as 

P{ Combat duration > t) = ( e -a ^ ) ( e~ B ^ ) . 

Now the probability that combat lasts until t and 
is terminated by a defender shooting down the 
bomber is 

P{Combat ends is ( dt ) , Defender wins) 

= (e _at adt )e _9t 

= (e-(° + B > t (a + #))(«/(»+•). 

This shows that a single combat duration is 
exp(a + b) and the event of a defender’s winning 
is independently a/(o + s). Likewise, the combat 
duration is exp(« + b) and a bomber’s win is, 
independently, b/ ( + b). If there are c combats 
going on then the first combat to ends does so 
in time exp(c(a + b)). 

Hence 

(b, c, d) ♦ (b, c-1 , d+1 ) with prob acdt 
(b, c, d) ■> (b + 1,c-1, d) with prob cdt . 
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(3) Simulation and Sojourn. The above shows that we 



may simulate the combat as follows. 

(i) You are in state (b, c, d). Obtain a 
sojourn time in that state that is 

exp(ibd + (a + e)c) 

i.e. ®bcd = 1 / (*bd + (a + 0 ) 0 ). 

The system stays in state (b,c,d) for 
time C®»S b Q d ) • 

(ii) With probability 

ibd / (ibd + (a + 0 ) 0 ) ♦ ( b — 1 , c + 1 , d-1) 

(NEW COMBAT) at time S bcd . 

(iii) With probability 

ac / (ibd + (a + 0 ) 0 ) ♦ (b, c — 1 , d+1 ) 

(DEFENDER SHOOTS DOWN BOMBER) at 
time S bcd . 

(iv) With probability 

0 C / (ibd + (a + 0 ) 0 ) ♦ (b + 1 ,c-1 ,d) 

(BOMBER SHOOTS DOWN DEFENDER) at 
time S bcd * 
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APPENDIX B 



FORTRAN PROGRAM LISTING FOR THE CRUDE SIMULATION 

OF THE BCD MODEL 



DIMENSION BX( 101 ) ,DX( 1 01 ) ,SEED(2) , BOX(2) ,D0X(2) 

INTEGER I , N , BB , DD , R 

REAL* 4 X , Y , Z , BX , DX , TXT , BOX , DOX 

DOUBLE PRECISION SEED 

DATA SEED / 1 234 . 0D0 , 567890 1 23 . 0D0/ 

R= 2 
N = 1 00 



RECEIVE INPUT DATA FROM TERMINAL 



8 



WRITE( * , 3 ) 

FORMAT ( 1 X , ’ENTER THE NUMBER OF BOMBERS’) 

READ ’ ( 12 ) ’ , BB 
WRITE( * ,4) 

FORMAT( 1 X , ’ENTER THE RATE WHICH A BOMBER SHOOTS 
DOWN A DEFENDER ’ ) 

READ ’ (F5 . 3 ) ’ , Y 
WRITE( * ,5) 

FORMAT ( 1 X , ’ENTER THE NUMBER OF DEFENDERS’) 

READ ’ ( 12 ) ’ , DD 
WRITE( * , 6 ) 

FORMAT( 1 X , ’ENTER THE RATE WHICH A DEFENDER’ 

& ’SHOOTS DOWN A BOMBER’) 

READ ’ (F5 . 3 ) ’ , X 

WRITE ( * , 7 ) ’ ENTER THE RATE WHICH FREE DEFENDERS’ 

& ’FIND FREE BOMBERS’ 

FORMAT ( 1 X , A) 

READ ’ (F5 • 3 ) ’ , Z 

WRITE( * ,8) ’ENTER THE TIME DURATION OF THE 
& ’BATTLE’ 

FORMAT ( 1 X , A) 

READ ’ (F5.3) ’ , TXT 



RUN REPLICATIONS OF N BATTLES AND OBTAIN SUMMARY OF 
N BATTLES 



CALL BATTLE( N , R , SEED , X , Y , Z , BB , DD , TXT , BX , DX ) 

COMPUTE STATISTICAL OUTPUT ANALYSIS OBTAIN PARAMETER 
ESTIMATES 



C 



CALL ST AT ( N , R , BX , DX , BOX , DOX ) 
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FORMAT AND PRINT OUTPUT OF PARAMETER ESTIMATES 
WRITE( 3,279) N 

279 FORMAT ( 15X , ’SAMPLE SI ZE ’ , 1 6 , / 1 5X , 6 ( ’ - ’ ) , IX , 

& 4( ) ) 

WRITE (3,280) BB , DD 

280 F0RMAT( //IX , 14 ,2X , ’BOMBERS’ ,6X , ’VERSUS’ , 6X , 
& I 4 , 2X , ’DEFENDERS’ ,/lX,13( ’-’ ),18X,14( ’-’ ) ) 

290 WRITE (3,300 )TXT 

300 F0RMAT( / / 18X , ’TIME’ ,F6.1,/,47(’-’)) 

WRITE( 3,310 ) 

310 FORMAT ( 1 9X , ’BOMBER’ ,2X, ’DEFENDER’ ,/18X, 

& 7 ( ’ = * ) ,2X,8( ’ = ’ ) ) 

WR I TE (3,320 ) BOX ( 1 ) , DOX ( 1 ) 

320 FORMAT (IX, ’AVERAGE’ ,7X,2F10.4) 

WRITE( 3 , 330 ) BOX ( 2 ) , DOX ( 2 ) 

330 FORMAT ( IX , ’ VARIANCE’ ,6X ,2F10 . 4 ) 

STOP 
END 



SUBROUTINE BATTLE 



SUBROUTINE 

& BATTLE ( N , R , SEED , X , Y , Z , BB , DD , TXT , BX , DX ) 
INTEGER BB , DD , I ,K ,N ,R 

REAL X ( N+l ) ,DX(N+1 ) , GA , SO J , X , Y , Z , B , C , D , NC , 
& BW,DW, INF, T, TIME, TXT 

DOUBLE PRECISION SEED(R) 

INITIALIZE STATISTICAL COUNTERS 

BX ( N+l )= 0.0 
DX ( N+l )= 0.0 

RUN N REPLICATIONS 



DO 200 I«1,N 

INITIALIZE START-TO-BATTLE VALUES 



T-0 . 0 
C =• 0 . 0 

B=REAL( BB ) 

D=REAL( DD ) 

BW- 0.0 
DW= 0.0 
NC= 1.0 

OBSERVE OCCURRENCE OF AN EVENT 
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1 00 



CALL UN IFOR ( SEED( \ ) , GA , 1 ) 



DETERMINE 

(B,C,D) 



NEXT INTERIM EVENT AND UPDATE THE STATE 



IF( GA .LE. BW) THEN 



B = B 


+ 


i . 0 


C = C 
D = D 




i . 0 


ELSE 


IF 


(GA .LE 


B = B 


- 


1 . 0 


C = C 


+ 


1 . 0 


D = D 
ELSE 
B = B 




1 . 0 


C = C 


- 


1 . 0 


D = D 


+ 


1 . 0 



LE. (BW+NC)) THEN 



END IF 

COMPUTE MEAN TIME IN STATE (B,C,D) 

IF ((B .EQ. 0.0 .OR. D .EQ. 0.0) .AND. C 
& .EQ. 0.0) THEN 

INF= 1 000000.0 

ELSE 

INF= 1 .0/(Z*D*B + (X + Y ) * C ) 

END IF 

GENERATE SOJOURN TIME IN STATE (B,C,D) 

CALL EXPON ( SEED( 2 ) , SO J , 1 ) 

T IME= -INF * ALOG(SOJ) 

ADVANCE THE SIMULATED TIME OF THE AIR BATTLE 
T= T + TIME 

COMPUTE PROBABILITY OF NEXT INTERIM EVENTS OCCURRING 



BW= 
DW= 
NC = 

CHECK CONDITIONS 



Y*C* INF 
X * C * INF 
Z*B*D* INF 

FOR OCCURRENCE 



OF TERMINATION EVENT 



IF (T .LT. TXT) GOTO 
ACCUMULATE SUMMARY OF N BATTLE 



1 00 

OUTCOMES 
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BX ( I )= B + C 
DX( I )= D + C 
BX ( N+ 1 )= BX ( N+ 1 ) + BX ( I ) 
DX (N + 1 ) = DX ( N+ 1 ) + DX ( I ) 
200 CONTINUE 
RETURN 
END 



SUBROUTINE EXPON 

SUBROUTINE EX PON ( SEED2 , A2 , K ) 

INTEGER I , K 

REAL ST IM , A2 , INF 

DOUBLE PRECISION EFF , SEED2 

EFF= 21 47483647 . 0D0 

SEED2=DM0D( 1 6807 . 0D0 * SEED2 , EFF ) 

A2 = SEED2/EFF 

RETURN 

END 



SUBROUTINE UNIFOR 

SUBROUTINE UNIFOR( SEED1 , A1 , K ) 
INTEGER I , K 
REAL A 1 

DOUBLE PRECISION EFF , SEED 1 

EFF= 21 47483647 . 0D0 

SEED1 =DMOD( 1 6807 . 0D0 * SEED1 , EFF ) 

A1 = SEED 1 /EFF 

RETURN 

END 



SUBROUTINE STAT 

SUBROUTINE STAT ( N , R , BX , DX , BOX , DOX ) 

INTEGER J , R , N 

REAL BX ( N + 1 ) ,DX(N+1 ) , BOX(R) ,DOX(R) 

BOX ( 2 ) = 0.0 
DOX ( 2 ) = 0.0 

COMPUTE THE ESTIMATES OF THE SAMPLE MEAN AND 
VARIANCE 

BOX ( 1 ) = BX (N + 1 )/N 
DOX ( 1 ) = DX ( N+1 )/N 
DO 260 1 = 1 ,N 

BOX ( 2 ) = BOX ( 2 ) + (BX(I)-BOX(I ) ) * * 2 
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260 



DOX ( 2 ) = DOX( 2 ) + (DX(- I )-DOX( 1 ) )**2 
CONTINUE 

BOX ( 2 ) = BOX (2)/ (N*( N- 1 ) ) 

DOX ( 2 ) = D0X(2)/(N*(N-1 )) 

RETURN 

END 
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APPENDIX C" 



FORTRAN PROGRAM LISTING OF THE BCD SIMULATION 
USING ANTITHETIC VARIATES 



DIMENSION BX( 1 01 ) , DX( 1 01 ) ,SEED( 2 ) ,BOX( 2 ) ,DOX( 2 ) 

INTEGER I,N,BB,DD,R 

REAL * 4 X , Y , Z , BX ,DX ,TXT , BOX , DOX 

DOUBLE PRECISION SEED 

DATA SEED / 1 234 . 0D0 , 567890 1 23 . 0D0/ 

R= 2 
N = 50 



RECEIVE INPUT DATA FROM TERMINAL 
WRITE( * ,3) 

FORMAT ( 1 X , ’ENTER THE NUMBER OF BOMBERS’) 

READ ’ ( 12 ) ’ , BB 
WRITE( * , 4 ) 

FORMAT( 1 X , ’ENTER THE RATE WHICH A BOMBER SHOOTS’ 
& ’DOWN A DEFENDER’) 

READ ’ (F5 • 3 ) ’ , Y 
WRITE( * ,5) 

5 FORMAT ( IX, ’ENTER THE NUMBER OF DEFENDERS’) 

READ ’ ( 12 ) ’ , DD 

WRITE( * , 6 ) 

6 FORMAT ( IX, ’ENTER THE RATE WHICH A DEFENDER 
& SHOOTS DOWN A BOMBER’) 

READ ’ (F5 . 3 ) ’ , X 

WRITE( * , 7 ) ’ ENTER THE RATE WHICH FREE DEFENDERS’ 
& ’FIND FREE BOMBERS’ 

7 FORMAT (IX, A) 

READ ’ (F5.3) ’ , Z 

WRITE( * , 8 ) ’ ENTER THE TIME DURATION OF THE’ 

& ’BATTLE’ 

8 FORMAT (IX, A) 

READ ’ ( F5 . 3 ) ’ , TXT 

RUN REPLICATIONS OF N BATTLES AND OBTAIN SUMMARY OF 
N BATTLES 

CALL BATTLE( N , R , SEED , X , Y , Z , BB , DD , TXT , BX , DX ) 

COMPUTE STATISTICAL OUTPUT ANALYSIS OBTAIN PARAMETER 
ESTIMATES 

CALL ST AT ( N , R , BX , DX , BOX , DOX ) 
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FORMAT AND PRINT OUTPUT OF PARAMETER ESTIMATES 



279 



280 

290 

300 

310 



320 

330 



WRITE( 3 , 279 ) N 

FORMAT( 1 5X , * SAMPLE SIZE’ , 16 , /I 5X , 6( ’ - ’ ) , 

& 1 X , 4 ( ’ - ’ ) ) 

WRITE( 3 , 280 ) BB.DD 

FORMAT ( / / 1 X , 14 , 2X , ’ BOMBERS ’ , 6X , ’ VERSUS ’ , 

& 6X , 14 , 2X , ’ DEFENDERS ’ ,/1X,13( ),1 8X , 1 4 ( » - * ) ) 

WRITE ( 3 , 300 )TXT 

FORMAT ( / / 1 8X , ’TIME’ , F6 . 1 , / , 47 ( ’ - ’ ) ) 

WRITE( 3,310) 

FORMAT ( 1 9X , ’ BOMBER * , 2X , ’ DEFENDER ’ , 

& / 1 8X , 7 ( ’ = ’ ),2X,8( ’ = ’ )) 

WRITE( 3,320) BOX ( 1 ) , DOX ( 1 ) 

FORMAT (IX,’ AVERAGE ’ , 7X , 2F 1 0 . 4 ) 

WRITE( 3,330) BOX ( 2 ) , DOX ( 2 ) 

FORMAT ( IX, ’VARIANCE’ ,6X,2F10.4) 

STOP 

END 



SUBROUTINE BATTLE 
SUBROUTINE 

& BATTLE( N , R , SEED , X , Y , Z , BB , DD , TXT , BX , DX ) 
INTEGER BB.DD, H, I, J,W,K,N,R 
REAL GA( 2 , 1 000 ) ,SOJ( 2, 1 000 ) , BX(N+1 ) , 

& DX ( N+ 1 ) , BAT( 50,2), DAT( 50,2), 

& X ,Y,Z,T, TXT, TIME, INF , NC , BW , DW , B , C , D 
DOUBLE PRECISION SEED(R) 

K= 1000 

INITIALIZE STATISTICAL COUNTERS 

BX ( N+ 1 )= 0.0 
DX ( N+ 1 )= 0.0 

RUN N REPLICATIONS 

DO 200 1=1 ,N 

CALL SOJOUR( SEED( 1 ) , SOJ , R , K ) 

CALL STATE( SEED( 2 ) , GA , R , K ) 

DO 175 J=1 , R 

INITIALIZE START-TO-BATTLE VALUES 



H = 0 
T = 0 . 0 
C = 0 . 0 

B=REAL( BB ) 
D=REAL ( DD ) 
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BW= 0.0 
DW= 0.0 
NC = 1 . 0 

OBSERVE OCCURRENCE OF NEXT EVENT 
00 H=H+1 

DETERMINE NEXT INTERIM EVENT AND UPDATE THE STATE 
( B , C , D ) 

IF(GA(J,H) .LE. BW) THEN 
B=B +1.0 
C=C - 1.0 
D = D 

ELSE IF (GA(J,H) . LE . ( BW+NC ) ) THEN 

B=B - 1.0 
C=C +1.0 
D=D - 1.0 
ELSE 
B = B 

C=C - 1.0 
D=D +1.0 
END IF 

COMPUTE MEAN TIME IN STATE (B,C,D) 

IF ((B .EQ. 0.0 .OR. D .EQ. 0.0) .AND. C 
& .EQ. 0.0) THEN 

INF= 1000000.0 

ELSE 

INF= 1 . 0 / ( Z *D * B + (X + Y ) *C ) 

END IF 

COMPUTE SOJOURN TIME IN STATE (B,C,D) 

TIME= -INF * ALOG( SO J ( J , H ) ) 

ADVANCE THE SIMULATED TIME OF THE AIR BATTLE 
T= T + TIME 

COMPUTE PROBABILITY OF NEXT INTERIM EVENTS OCCURRING 



BW= Y*C*INF 
DW= X * C * INF 
NC= Z*B*D*INF 

CHECK FOR OCCURRENCE OF TERMINATION EVENT 
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IF (T .LT. TXT) GOTO 100 
C 

C RECORD RESULTS OF BATTLE 
C 

BAT ( I , J ) = B + C 
DAT ( I , J ) = D + C 
175 CONTINUE 

C 

C ACCUMULATE SUMMARY OF N BATTLE OUTCOMES 
C 

BX ( I ) = ( BAT (1,1) + BAT (1,2)) *.5 
DX ( I ) = ( DAT (1,1) + DAT ( I , 2 ) ) * . 5 
BX ( N+ 1 ) = BX ( N+ 1 ) + BX(I,J) 

DX ( N+ 1 ) = DX ( N+1 ) + DX ( I , J ) 

200 CONTINUE 
RETURN 
END 
C 
C 

C SUBROUTINE SOJOUR 
C 

SUBROUTINE SO JOUR ( SEED2 , A2 , W , K ) 
INTEGER I,W,K 
REAL A2(W,K) 

DOUBLE PRECISION EFF , SEED2 
EFF= 21 47483647. 0D0 
DO 10 1=1 , K 

SEED2=DMOD( 1 6807 . 0D0 * SEED2 , EFF ) 
A2( 1 , I ) = SEED2/EFF 
A2 ( 2 , I ) = 1 .0 - A2( 1 ,1) 

10 CONTINUE 
RETURN 
END 
C 
C 

C SUBROUTINE STATE 
C 

SUBROUTINE STATE ( SEED1 ,A1 ,W,K) 
INTEGER I,K,W 
REAL A 1 ( W , K ) 

DOUBLE PRECISION EFF , SEED1 
EFF= 21 47483647 . 0D0 
DO 10 1 = 1 ,K 

SEED1 =DMOD( 1 6807 . 0D0 * SEED 1, EFF) 
A1 ( 1 , I )= SEED1 /EFF 
A1 ( 2 , I ) = 1.0 - A1 ( 1 , I ) 

10 CONTINUE 
RETURN 
END 
C 
C 
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SUBROUTINE STAT 

SUBROUT INE STAT ( N , R , BX , DX , BOX , DOX ) 

INTEGER J , R , N 

REAL BX ( N+1 ) ,DX(N+1 ) , BOX(R) ,DOX(R) 

BOX ( 2 ) = 0.0 
DOX ( 2 ) = 0.0 

COMPUTE THE ESTIMATES OF THE SAMPLE MEAN AND 
VARIANCE 

BOX ( 1 ) = BX ( N+1 , J)/N 
DOX ( 1 ) = DX (N+1 , J)/N 
DO 260 1 = 1 ,N 

BOX ( 2 ) = BOX ( 2 ) + (BX(I)-BOX(I ))**2 
DOX ( 2 ) = DOX ( 2 ) + (DX( I )-DOX( 1 ) )**2 
260 CONTINUE 

BOX ( 2 ) = B0X(2)/(N*(N-1 )) 

DOX ( 2 ) = D0X(2)/(N*(N-1 )) 

RETURN 

END 
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APPENDIX D 



FORTRAN PROGRAM LISTING OF THE BCD SIMULATION USING 

STRATIFIED SAMPLING 



DIMENSION BX( 1 01 ) ,DX( 1 01 ) ,SEED(2 ) , BOX( 2 ) ,DOX(2 ) 

INTEGER I , N , BB , DD , R 

REAL*4 X,Y,Z,BX,DX,TXT, BOX , DOX 

DOUBLE PRECISION SEED 

DATA SEED / 1 234 . 0D0 , 567890 1 23 . 0D0/ 

R= 2 
N=1 00 



RECEIVE INPUT DATA FROM TERMINAL 
WRITE ( * ,3) 

FORMAT ( 1 X , ’ENTER THE NUMBER OF BOMBERS’) 

READ ’ ( 12 ) ’ , BB 
WRITE( * , 4 ) 

FORMAT ( 1 X , ’ENTER THE RATE WHICH A BOMBER SHOOTS’ 
& ’DOWN A DEFENDER’) 

READ ’ ( F5 . 3 ) ’ , Y 
WRITE( * , 5 ) 

FORMAT ( 1 X , ’ENTER THE NUMBER OF DEFENDERS’) 

READ ’ ( 12 ) ’ , DD 
WRITE( * , 6 ) 

FORMAT ( 1 X , ’ENTER THE RATE WHICH A DEFENDER’ 

& ’SHOOTS DOWN A BOMBER’) 

READ ’ (F5.3 ) ’ , X 

WRITE( * , 7 ) ’ ENTER THE RATE WHICH FREE DEFENDERS’ 
& ’FIND FREE BOMBERS’ 

7 FORMAT ( 1 X , A ) 

READ ’ (F5.3) ’ , Z 

WR ITE( * ,8) ’ ENTER THE TIME DURATION OF THE’ 

& ’BATTLE’ 

FORMAT ( 1 X ,A) 

READ ’ (F5 • 3 ) ’ , TXT 

RUN REPLICATIONS OF N BATTLES AND OBTAIN SUMMARY OF 
N BATTLES 

CALL BATTLE( N , R , SEED , X , Y , Z , BB , DD , TXT , BX , DX ) 

COMPUTE STATISTICAL OUTPUT ANALYSIS OBTAIN PARAMETER 
ESTIMATES 

CALL STAT ( N , R , BX , DX , BOX , DOX ) 
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c 

C FORMAT AND PRINT OUTPUT OF PARAMETER ESTIMATES 
C 

279 



280 



290 

300 

310 



320 

330 



C 
C 

C SUBROUTINE BATTLE 
C 

SUBROUTINE 

& BATTLE( N , R , SEED , X , Y , Z , BB , DD , TXT , BX , DX ) 
INTEGER BB , DD ,H,I,J,K,N,R,W 
REAL GA( 3 , 1 000 ) ,SOJ( 3 ,1000) ,BM( 34,3), 

& DF( 34 , 3 ) , BX ( N+1 ),DX(N+1 ), X , Y , Z , T , TXT , TIME , 
& INF , NC , BW , DW , B , C , D , BM , DF 
DOUBLE PRECISION SEED(R) 

K= 1000 
W= 3 
C 

C INITIALIZE STATISTICAL COUNTERS 
C 

BX (N+1 )= 0.0 
DX (N+1 ) = 0.0 
C 

C RUN N REPLICATIONS 
C 

DO 200 1 = 1 , N 

CALL SO J OUR( SEED( 1 ) , SO J , W , K ) 

CALL STATE( SEED( 2 ) , GA , W , K ) 

DO 175 J = 1 ,W 
C 

C INITIALIZE START-TO-BATTLE VALUES 
C 

H = 0 
T = 0 . 0 



WRITE( 3 , 279 ) N 

FORMAT ( 1 5X , ’ SAMPLE SIZE’ , I 6 , / 1 5X , 6 ( ), 

& 1 X , 4 ( ’ - ’ ) ) 

WRITE( 3 , 280 ) BB , DD 

& FORMAT (//IX, 14 , 2X , ’BOMBERS’ ,6X, ’VERSUS’ , 
& 6X, 14, 2X, ’DEFENDERS’ , 

& /IX, 13( ) , 18X, 1 4 ( ’-’ ) ) 

WRITE (3,300)TXT 

FORMAT ( / / 1 8X , ’TIME’ ,F6. 1 ,/,47( ’-’ )) 

WR ITE( 3,310) 

FORMAT ( 19X, ’BOMBER’ , 2X , ’DEFENDER’ , 

& / 1 8X , 7 ( ’ = ’ ) ,2X ,8( ’ = ’ ) ) 

WRITE( 3 , 320 ) BOX ( 1 ) , DOX ( 1 ) 

FORMAT ( IX, ’AVERAGE’ ,7X,2F10.4) 

WRITE( 3,330) BOX ( 2 ) , DOX ( 2 ) 

FORMAT ( IX, ’VARIANCE’ ,6X,2F10.4) 

STOP 

END 
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C = 0.0 

B=REAL( BB ) 

D=REAL( DD ) 

NC= 1 . 0 
BW= 0 . 0 
DW= 0 . 0 

OBSERVE AN OCCURRENCE OF AN EVENT 
100 H=H+1 

DETERMINE NEXT INTERIM EVENT AND UPDATE THE STATE 
( B , C , D ) 

IF ( GA( J , H ) .LE. BW) THEN 
B=B +1.0 
C=C - 1.0 
D = D 

ELSE IF ( GA ( J , H ) . LE . ( BW+NC ) ) THEN 

B=B - 1.0 
C=C +1.0 
D=D - 1.0 
ELSE 
B = B 

C=C - 1.0 
D=D+ 1 . 0 
END IF 

COMPUTE MEAN TIME IN STATE (B,C,D) 

IF ( ( B .EQ. 0.0 .OR. D .EQ. 0 . 0 ) . AND . C 
& .EQ. 0.0) THEN 

INF= 1000000.0 
ELSE 

INF= 1 .0/(Z*B*D + (X+Y)*C) 

END IF 

COMPUTE SOJOURN TIME OF STATE (B,C,D) 

TIME = -INF * ALOG( SO J ( J , H ) ) 

ADVANCE SIMULATED TIME OF THE AIR BATTLE 
T= T + TIME 

COMPUTE PROBABILITY OF NEXT INTERIM EVENTS OCCURRING 
C 

NC= Z*B*D*INF 
BW= Y * C * INF 
DW= X*C* INF 
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CHECK FOR OCCURRENCE OF TERMINATION EVENT 

IF (T .LT. TXT ( TI ) ) GOTO 100 

RECORD RESULTS OF BATTLE 

BM( I , J ) = B + C 
DF( I , J ) = D + C 
175 CONTINUE 

ACCUMULATE SUMMARY OF N BATTLE OUTCOMES 

BX ( I ) = ( BM( 1,1) + BM( 1,2) + BM( I , 3 ) / 3 • 0 
DX ( I ) = ( DF (1,1) + DF (1,2) + BM(I,3)/3.0 
BX(N+1 )« BX ( N+ 1 ) + BX ( I ) 

DX ( N+1 ) = DX (N+1 ) + DX ( I ) 

200 CONTINUE 
RETURN 
END 



SUBROUTINE SOJOUR 

SUBROUTINE SO JOUR ( SEED2 , B2 , W , K ) 

INTEGER I , W , K , J 

REAL B2 ( W , K ) , A2 

DOUBLE PRECISION EFF , SEED2 

EFF= 21 47483647 . 0D0 

DO 10 1=1 ,K 

SEED2=DM0D( 1 6807 . 0D0 * SEED2 , EFF ) 

A2= SEED2/EFF 
DO 5 J=1 ,W 

B2 ( J , I ) = AM0D(A2 + ( ( J - 1 ) * 1 . 0 ) / 3 . 0 , 1 . 0 ) 
5 CONTINUE 

10 CONTINUE 
RETURN 
END 



SUBROUTINE STATE 

SUBROUTINE STATE( SEED 1 , A1 , W , K ) 

INTEGER I , K , W , J 

REAL A1 ( W , K ) , A2 

DOUBLE PRECISION EFF , SEED1 

EFF= 2147483647. 0D0 

DO 10 1 = 1 ,K 

SEED1 =DMOD( 1 6807 . 0D0 * SEED1 , EFF ) 
A2= SEED 1 /EFF 
DO 5 J= 1 ,W 
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1 .0)/3.0, 1 .0) 



A 1 ( J , I ) = AM0D(A2 + ( (- J - 1 ) * 
5 CONTINUE 

RETURN 
END 



SUBROUTINE STAT 

SUBROUTINE STAT ( N , R , BX , DX , BOX , DOX ) 

INTEGER R , N 

REAL BX ( N+ 1 ) ,DX(N+1 ) , BOX(R) ,DOX(R) 

BOX ( 2 ) = 0.0 
DOX ( 2 ) = 0.0 

COMPUTE THE ESTIMATES OF THE SAMPLE MEAN AND 
VARIANCE 

BOX ( 1 ) = BX ( N+ 1 )/N 
DOX ( 1 ) = DX( N+1 )/N 
DO 260 1 = 1 ,N 

BOX ( 2 ) = BOX ( 2 ) + (BX(I)-BOX(I ) ) * * 2 
DOX ( 2 ) = DOX ( 2 ) + (DX( I )-DOX( 1 ) )**2 
260 CONTINUE 

BOX ( 2 ) = BOX (2)/(N*(N-1 )) 

DOX ( 2 ) = DOX (2)/(N*(N-1 )) 

RETURN 

END 
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APPENDIX E 



STATISTICAL OUTPUT DATA FROM THE BCD SIMULATIONS 



TABLE E . 1 OUTPUT PARAMETERS OF THE BCD MODEL ESTIMATED 
FROM CRUDE SIMULATION 



SIMULATION 


SCENARIO 




DEFENDER 


BOMBER 


DEFENDER 


BOMBER 


TIME 


MEAN 


VAR 


MEAN 


VAR 


1 


1 0 


1 0 


25 


7 . 6 


.0191 


7.5 


. 0250 


2 


1 0 


1 0 


75 


4.4 


. 0387 


4.5 


. 0358 


3 


1 0 


1 0 


1 25 


2.8 


.0375 


3.4 


.0417 


4 


1 0 


30 


25 


5.3 


.0328 


25.1 


. 0541 


5 


1 0 


30 


75 


1 .3 


.0157 


21.3 


. 1 402 


6 


1 0 


30 


1 25 


.3 


. 0036 


20.8 


. 1 725 


7 


1 0 


50 


25 


4.2 


. 0251 


44.7 


. 0558 


8 


1 0 


50 


75 


. 6 


. 0064 


40.8 


. 1 654 


9 


1 0 


50 


1 25 


. 1 


.0014 


39.7 


. 1 874 


1 0 


30 


1 0 


25 


25 . 1 


. 0458 


5.4 


. 0350 


1 1 


30 


1 0 


75 


21 . 1 


. 1 407 


1 .4 


.0141 


1 2 


30 


1 0 


1 25 


19.8 


.2247 


. 4 


. 0045 


1 3 


30 


30 


25 


18.7 


.0933 


18.5 


.0842 


1 4 


30 


30 


75 


8.2 


.1618 


8.4 


.1895 


1 5 


30 


30 


1 25 


5.7 


. 1 984 


4.7 


.1341 


1 6 


30 


50 


25 


14.4 


. 0628 


34.5 


.1104 


1 7 


30 


50 


75 


3.4 


.0362 


22 . 9 


.3479 


1 8 


30 


50 


1 25 


.9 


.0136 


20 . 1 


. 4881 


1 9 


50 


1 0 


25 


44.4 


. 0651 


4 . 1 


.0317 


20 


50 


1 0 


75 


40.7 


. 1 427 


.7 


. 0063 


21 


50 


1 0 


1 25 


41 . 1 


. 1 666 


. 1 


.0014 


22 


50 


30 


25 


34.3 


. 1 687 


14.0 


. 0826 


23 


50 


30 


75 


23.5 


.3339 


3. 1 


. 0359 


24 


50 


30 


1 25 


21 . 1 


. 5300 


.8 


.0107 


25 


50 


50 


25 


27.9 


. 1 590 


26.4 


.1501 


26 


50 


50 


75 


9.5 


.2717 


11.5 


.2865 


27 


50 


50 


1 25 


5.2 


.2126 


7.9 


.2986 
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TABLE E . 2 OUTPUT PARAMETERS OF THE BCD MODEL ESTIMATED 
FROM THE SIMULATION USING ANTITHETIC VARIATES 



SIMULATION 


SCENARIO 
DEFENDER BOMBER 


TIME 


DEFENDER 
MEAN VAR 


BOMBER 
MEAN VAR 


1 


10 


10 


25 


7.6 


. 0067 


7 . 5 


. 0105 


2 


10 


10 


75 


4 . 5 


. 0048 


4 . 4 


. 0060 


3 


10 


10 


125 


3 . 2 


. 0079 


3 . 3 


.0057 


4 


10 


30 


25 


5 . 1 


. 0129 


25 . 1 


. 0100 


5 


10 


30 


75 


1 . 2 


. 0071 


21 . 3 


. 0329 


6 


10 


30 


125 


. 3 


. 0024 


20 . 2 


. 0446 


7 


10 


50 


25 


4 . 4 


.0123 


44 . 4 


.0146 


8 


10 


50 


75 


. 6 


. 0060 


40 . 9 


. 0395 


9 


10 


50 


125 


. 1 


.0013 


40 . 1 


. 0683 


10 


30 


10 


25 


25 . 4 


. 0116 


5 . 0 


. 0116 


11 


30 


10 


75 


21 . 1 


. 0271 


1 . 2 


. 0056 


12 


30 


10 


125 


20.2 


. 0558 


. 3 


.0030 


13 


30 


30 


25 


18.5 


. 0164 


18.4 


. 0210 


14 


30 


30 


75 


8.5 


. 0251 


8.4 


. 0270 


15 


30 


30 


125 


5 . 1 


. 0307 


5 . 1 


. 0352 


16 


30 


50 


25 


14 . 6 


. 0223 


34 . 6 


. 0308 


17 


30 


50 


75 


3.0 


. 0138 


23 . 6 


. 0542 


18 


30 


50 


125 


1 . 0 


. 0075 


20 . 6 


. 0446 


19 


50 


10 


25 


44 . 0 


.0169 


4 . 3 


.0104 


20 


50 


10 


75 


40 . 2 


. 0410 


. 6 


. 0059 


21 


50 


10 


125 


40 . 6 


. 0451 


. 1 


.0015 


22 


50 


30 


25 


34 . 3 


. 0279 


14 . 4 


. 0368 


23 


50 


30 


75 


23 . 1 


. 0375 


3 . 1 


. 0194 


24 


50 


30 


125 


20 . 5 


. 0739 


. 9 


. 0077 


25 


50 


50 


25 


27 . 3 


. 0266 


27.1 


. 0357 


26 


50 


50 


75 


10 . 5 


. 0269 


10 . 7 


. 0309 


27 


50 


50 


125 


6 . 4 


. 0458 


6 . 4 


. 0472 



102 



1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 

21 

22 

23 

24 

25 



EFFICIENCY OF ANTITHETIC VARIATES IN THE BCD 



SCENARIO DEFENDER BOMBER 



DEFENDER 


BOMBER 


TIME 


RE 


IP(*) 


RE 


IP(£) 


10 


10 


25 


2 . 9 


65 . 0 


2 . 4 


58 . 0 


10 


10 


75 


8.1 


87 . 6 


6.0 


83.2 


10 


10 


125 


4.7 


78.9 


7 . 3 


86.3 


10 


30 


25 


2 . 5 


60.7 


5 . 4 


81 . 5 


10 


30 


75 


2 . 2 


54 .8 


4 . 3 


76 . 5 


10 


30 


125 


1 . 5 


33.3 


3.9 


74 . 1 


10 


50 


25 


2 . 0 


51 . 0 


3.8 


73.8 


10 


50 


75 


1 . 1 


6.3 


4 . 2 


76 . 1 


10 


50 


125 


1 . 1 


7 . 1 


2.7 


63.6 


30 


10 


25 


3.9 


74 . 7 


3 . 0 


66 . 9 


30 


10 


75 


5 . 2 


80.7 


2 . 5 


60.3 


30 


10 


125 


4 . 0 


75 . 2 


1 . 5 


33.3 


30 


30 


25 


5 . 7 


82 . 4 


4 . 0 


75 . 1 


30 


30 


75 


6 . 4 


84 . 5 


7 . 0 


85.8 


30 


30 


125 


6 . 5 


84.5 


3.8 


73.8 


30 


50 


25 


2.8 


64 . 5 


3.6 


72 . 1 


30 


50 


75 


2 . 6 


61 . 9 


6 . 4 


84 . 4 


30 


50 


125 


1.8 


44 . 9 


10 . 9 


90 . 9 


50 


10 


25 


3.9 


74 . 0 


3 . 0 


67 . 2 


50 


10 


75 


3.5 


71 . 2 


1 . 1 


6.3 


50 


10 


125 


3.7 


72 . 9 


. 9 


7 . 1 


50 


30 


25 


6 . 0 


83.4 


2 . 2 


55.4 


50 


30 


75 


8.9 


88.8 


1 . 9 


46 . 0 


50 


30 


125 


7 . 2 


86 . 1 


1 . 4 


28.0 


50 


50 


25 


6.0 


83.3 


4 .2 


76.2 


50 


50 


75 


10 . 1 


90 . 1 


9.3 


89 . 2 


50 


50 


125 


4 . 6 


78.5 


6.3 


84 . 2 
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TABLE E . 4 OUTPUT PARAMETERS OF THE BCD MODEL ESTIMATED 
FROM SIMULATION USING STRATIFIED 'SAMPLING 



SIMULATION 


SCENARIO 




DEFENDER 


BOMBER 


DEFENDER 


BOMBER 


TIME 


MEAN 


VAR 


MEAN 


VAR 


1 


1 0 


1 0 


25 


7.4 


.0119 


7.5 


.0114 


2 


1 0 


1 0 


75 


4.4 


. 0248 


4.5 


. 0237 


3 


1 0 


1 0 


1 25 


3.1 


.0160 


3 . 2 


.0197 


4 


1 0 


30 


25 


5.1 


.0124 


25.0 


.0252 


5 


1 0 


30 


75 


1 . 2 


. 0076 


21.2 


. 0710 


6 


1 0 


30 


1 25 


.4 


. 0038 


20 . 5 


. 0801 


7 


1 0 


50 


25 


4.4 


.0126 


44 . 0 


.0375 


8 


1 0 


50 


75 


.9 


. 0083 


40.2 


.1021 


9 


1 0 


50 


1 25 


. 1 


.0006 


40 . 0 


. 0930 


1 0 


30 


1 0 


25 


25 . 1 


.0233 


5.3 


.0124 


1 1 


30 


1 0 


75 


21.1 


. 0726 


1 . 3 


.0118 


1 2 


30 


1 0 


1 25 


20 . 3 


. 0877 


.4 


. 0027 


1 3 


30 


30 


25 


18.4 


. 0283 


18.4 


.0193 


1 4 


30 


30 


75 


8 . 0 


. 0439 


8 . 2 


. 0700 


15 


30 


30 


1 25 


5.5 


. 0505 


4.8 


. 0426 


1 6 


30 


50 


25 


14.2 


.0219 


34.5 


. 0798 


1 7 


30 


50 


75 


3.2 


.0211 


23.3 


. 1 236 


1 8 


30 


50 


1 25 


.8 


. 0081 


20 . 1 


. 1 880 


1 9 


50 


1 0 


25 


44.6 


.0339 


4.3 


.0163 


20 


50 


1 0 


75 


41.0 


. 0872 


.7 


. 0053 


21 


50 


1 0 


1 25 


40 . 4 


. 1 089 


. 1 


.0012 


22 


50 


30 


25 


34.6 


. 0571 


14.3 


. 0439 


23 


50 


30 


75 


22.6 


.1187 


3.3 


.0197 


24 


50 


30 


1 25 


20 . 9 


. 2402 


1 . 0 


. 0086 


25 


50 


50 


25 


27.7 


. 0630 


27 . 0 


. 1 084 


26 


50 


50 


75 


10.5 


. 0898 


10.7 


.0881 


27 


50 


50 


1 25 


6.5 


.0911 


6.4 


. 1 084 
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1 

2 

3 

4 

5 

6 

7 

8 

9 

1 0 

1 1 

1 2 

13 

1 4 

1 5 

1 6 

1 7 

18 

1 9 

20 

21 

22 

23 

24 



EFFICIENCY OF STRATIFIED SAMPLING IN THE BCD 



SCENARIO 




DEFENDER 


BOMBER 


DEFENDER 


BOMBER 


TIME 


RE 


IPW 


RE 


IP(*) 


1 0 


1 0 


25 


1 . 6 


39 . 0 


2 . 2 


54 . 4 


1 0 


1 0 


75 


1 . 6 


36 . 0 


1 .5 


33.8 


1 0 


1 0 


1 25 


2.3 


57.3 


2 . 1 


52.8 


1 0 


30 


25 


2 . 6 


62 . 2 


2 . 1 


53.4 


1 0 


30 


75 


2 . 1 


51.6 


2 . 0 


49 . 4 


1 0 


30 


1 25 


.9 


-5.5 


2.2 


53 . 6 


1 0 


50 


25 


2 . 0 


49.8 


1 .5 


32.8 


1 0 


50 


75 


.8 


-29 . 7 


1 . 6 


38.3 


1 0 


50 


1 25 


2.3 


57 . 1 


2 . 0 


50 . 4 


30 


1 0 


25 


2 . 0 


49 . 1 


2.8 


64.6 


30 


1 0 


75 


1 .9 


48.4 


1 . 2 


16.3 


30 


1 0 


1 25 


2 . 6 


61.0 


1 . 7 


40 . 0 


30 


30 


25 


3.3 


69.7 


4.4 


77 . 1 


30 


30 


75 


3.7 


72 . 9 


2.7 


63 . 1 


30 


30 


1 25 


3.9 


74 . 5 


3 . 1 


68 . 2 


30 


50 


25 


2.9 


65 . 1 


1 .4 


27 . 7 


30 


50 


75 


1 .7 


41.7 


2.8 


64 . 5 


30 


50 


1 25 


1 . 7 


40 . 7 


2 . 6 


61.5 


50 


1 0 


25 


1 .9 


47.9 


1 .9 


48.6 


50 


1 0 


75 


1 . 6 


39.0 


1 . 2 


15.9 


50 


1 0 


125 


1 .5 


34 . 6 


1 . 2 


14.2 


50 


30 


25 


3.0 


66 . 2 


1 . 9 


46 . 9 


50 


30 


75 


2.8 


64.5 


1 .8 


45.1 


50 


30 


1 25 


2 . 2 


54.7 


1 . 2 


19.6 


50 


50 


25 


2.5 


60 . 4 


1 .4 


27.8 


50 


50 


75 


3.0 


66.9 


3.3 


69 . 2 


50 


50 


1 25 


2.3 


57. 1 


2.8 


63 . 7 
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